{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import utils\n",
    "\n",
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Sensitivity, Specificity, Precision"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Q1. Given a test result is positive, what's the probability that the patient is really positive?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From the Bayesian's rule:\n",
    "\n",
    "$$Pr(real=P|test=P) = \\frac{Pr(test=P|real=P)*Pr(real=P)}{Pr(test=P)} = \\frac{Pr(test=P|real=P)*Pr(real=P)}{Pr(test=P|real=P)*Pr(real=P) + Pr(test=P|real=N)*Pr(real=N)}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Therefore,\n",
    "$$ Pr(real=P|test=P) = \\frac{Sensitivity * Prior}{Sensitivity * Prior + (1-Specificity)*(1-Prior)} \\;\\;\\;\\;\\;\\;(1)$$|"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As of 4/21/2020, IL \"confirmed\" 33k, with population 12.67M, So the prior \n",
    "$$Pr(P)\\approx 33k/12.67M = 0.26\\% $$\n",
    "\n",
    "As of 4/24/2020, NY \"confirmed\" 271k, with population 19.45M, So the prior \n",
    "$$Pr(P)\\approx 271k/19.45M = 1.39\\% $$\n",
    "\n",
    "In the whole U.S., \"confirmed\" 827k, with population 328M, so the prior\n",
    "$$Pr(P)\\approx 827k/328M = 0.25\\% $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$Sensitivity = Pr(test=P | real=P)$$\n",
    "$$Specificity = Pr(test=N | real=N)$$\n",
    "$$Prior = Pr(real=P)$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$ Pr(real=P|test=P) \\approx \\frac{80.3\\%*.05\\%}{80.3\\% * .05\\% + (1-99.5\\%)*(1-.05\\%)} \\approx 8.5\\%$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "prior = 0.002"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def posterior_pos(sen, spe, prior):\n",
    "    return sen*prior / (sen * prior + (1-spe) * (1 - prior))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0005789473684210527"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "1.1/1900"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 202,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.679003512423572"
      ]
     },
     "execution_count": 202,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "posterior_pos(.803, .995, prior=.013)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = np.arange(1e-4, .3, 1e-4) * 100\n",
    "sens = 0.99\n",
    "frm = pd.DataFrame({'prior (%)': x, \n",
    "                    '99.8% specificity': posterior_pos(sens, .998, x/100) * 100,\n",
    "                    '99.5% specificity': posterior_pos(sens, .995, x/100) *100, \n",
    "                   '99% specificity': posterior_pos(sens, .99, x/100) * 100})\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x19494dece80>"
      ]
     },
     "execution_count": 53,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtQAAAG9CAYAAAA1PuH/AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOzdd3gU1frA8e9JsumVFCCEEAKhBpLQQSQoVUAQLKDgtdDVC3iVfuUiouLFclG8cMUfFqQJIqIgIiIC0nsLLRAgJKSS3ja75/fHTEICCdKSTTmf59lnN1PfmR3Yd868c0ZIKVEURVEURVEU5e5YWToARVEURVEURanMVEKtKIqiKIqiKPdAJdSKoiiKoiiKcg9UQq0oiqIoiqIo90Al1IqiKIqiKIpyD1RCrSiKoiiKoij3QCXUilKJCSGWCyEeK6NlbxVCjCiLZZewLn8hRIYQwvoW02QIIQLLYN1RQoiA+73cikwI8bMQ4rlbjF8ohHjjPqxnjRCi970upyIRQkwTQnx+i/FDhRCbbnNZJ4QQXe9bcIqiWIxKqBWlgtATu2w9cYwTQnwhhHC+xfQtgRDgh/KLsnDdM4UQRj3WFCHETiFEx7tdnpTykpTSWUpp0pd/UzKvjz9/r7HfCSHE80KIHfdpWVFCiO53MV+AEELq+zpDX86Ue4lFSvmIlPIrffk3baOUcoyU8q17WYduDvD2X00khHAUQkwSQhwQQiQLIS4LIZYKIcL+Yj5bIcQHQohofd9cEEJ8dB/iLpWU8h0p5Qh9/QXfjU2R8UullD1vc1nNpZRb9WXNFEJ8UyZB34IQoqkQYosQIlUIcU4IMfCG8SP04RlCiI1CCN8i454RQsTq+71rkeEN9P8TSj1BVpSqRiXUilKxPCqldAZaAW2Bf944QZEf79HAUlnK05mK/siXkZV6rN7ADmCNEEKU8TqrM3d9fz8NzKgMLb9Syr2AqxCiTWnTCCFqAbuABsCLQC2gKbAG+EYI8cItVjEVaAO0A1yAh4BD9yf6qk//P+IH4CegBjAKbZ830seHA+8AA/TxF4DlReadg/Z/1d+B+UUW/THwj4ITZEWpDlRCrSgVkJTyCvAzEAygt4K9LIQ4C5zVJ3sE+KNgHr2l8U8hxEdCiGRgphDCTgjxvhDikt7qvVAI4aBP7yGE+EkIkSCEuKZ/9ruLWI3AV2iJkKcQwlcIsU5vaTwnhBhZJMZ2Qoj9Qog0PZ4P9eGFLX1CiLeBB4H5eqvY/CL7oKEQooMQ4mrR1i8hxEAhxFH9s5UQYooQIlIIkSSE+FYIUeNOt0sI0RRYCHQsaInXh99qn3rp+zFF3/7tejxLAH/gR31Zk+40ngJSyl3ACSBYX/Y/hRAXhRDxQoivhRBueiz2Qohv9H2QIoTYJ4SoqY/bqrc8lraNXwohZuufI4QQ/YrsFxshRKIQopX+dwe9NTJFCHFE3FzCsBXoe4tNWgbMk1KOllIekVLmSSkzpJTfAZ2B14QQDUqZty3wvZQyRmqipJRfF4nVVwjxnX6MXxBCjCsybqZ+bHwthEgXWvlFmyLjJwshrujjTgshuhWZr6AleZv+nqLvv46iSIu/fmy8XzRgIcQPQoh/6J+jhBDdhXZyNA0YrC/niBDiSSHEgRvmfU0IsfYW+/JONQF8gY+klCYp5RbgT+BZffyjwCop5QkpZR7wFtBF/z48gStSylhgMxCox/iEPnz3fYxTUSo8lVArSgUkhKgL9KF4a9tjQHugmRDCCagPnL5h1vbAecAH7VL7e0AjIBRoCNQBZujTWgFfAPXQkr1sircy3W6sdsDzQLSUMhGtBSsa7Yf6CeCdgmQEmIeWPLmitUh+e+PypJTTge3AK3qZxys3jN8NZAIPFxn8DFpiBjAObV+F6zFcAz690+2SUkYAY4Bdehzu+qhb7dPX9G33BmqiJUlSSvkscAn9CoSU8t8AehJa2uumsg6heQBojnZsPK+/HkJLaJy5/h0+B7gBddGSnzFo3/HtbGNRy9FaxQv0AhKllAeFEHWA9cBstBbM14HvhBDeRaaPQCtNuonQWkCNUsrFQoi6Qis9iBFCfCqE2CulvAa8C4wtaX5gN/APIcRLQogWQly/QiKEsAJ+BI6gfUfdgAlCiF5F5u8PrADcgXXo+04I0Rh4BWgrpXTRtzmqhPV30d/d9f2364bxy9CSZKEv1wPoqa+zkJRyI1pL8Ep9OSF6PPX1k54Cw4AlJe0IIcR/b3EsHS1pHqCkK0oC/URe/yxuGIc+PgHtBNoP6AGcEFqJ2j/RrhwoSrWiEmpFqVjW6q2EO9Ban98pMu5dKWWylDIbLQEASL9h/hgp5SdSynwgBxgJvKrPl64vbwiAlDJJSvmdlDJLH/c2WhJ6u57SY70MtAYe008EOgOTpZQ5UsrDwOdcb/EyAg2FEF56K+TdtmIVJnlCCBe0k4/l+rjRwHQpZbSUMheYCTwh7kMJjJ4YlbpP0bavNlBPSmmUUm4vrSQHQErpfovXnBsmTwSS0fbnFCnlb8BQ4EMp5XkpZQZaIjNE31YjWiLdUG99PCClTLuLzV4G9BdCOOp/Fz15GQZskFJukFKapZS/AvvRvo8C6Vw/Xm/Ug+vJ5fvATrSTu5/QSjkADqO1pJbkXbQTnKH6eq+I6zdbtgW8pZSz9Fbv88Airn9XADv02E1oiWpB4m8C7NBOXg16y3dkKTHcynZAol1xAe0Ec5eUMuavZtSP3ZVo+xghRHMgAG3flDT9S7c4llqWsppTQDwwUQhhEEL0RPs/oOC73oD277yl0K7CzNC3x1FKaUY70VmNdiI1EpgFfAK0EEL8LoT4RQgRjKJUAyqhVpSK5TH9B7Ce/gNZtEXxcpHPKfq7yw3zF53GG+2H8UBBSxWwUR9ecCPY/4RWLpCGdvnaXdz+jUTf6rH6SCkfllIeQGsRLkg0C1xEayEEGI7WuntKaCUI/bg7y4BBeuv4IOCglPKiPq4e8H2RbY5AS5Bq3uW6irrlPgXmAueATUKI8yW1Mt8DLymlh5SyqZTyY32YL9r+LXARsEHb1iXAL8AKvdX330IIw52uVEp5Dm0fPqon1f25nlDXA54s2hqKdkJVu8giXLh+vN7IB7iif24BLJNS5kspf0Y7gQCthf1KSTPrJwqfSikfQEva3wYW66269QDfG2KbRvHj4GqRz1mAvRDCRt/mCWgnY/FCiBWiyM14t0s/mVrB9Rb+Z4Cld7CIr4Bn9BO5Z9H+zeXeaRy3iM+IdjWnL9q+eA3tqlG0Pv434F/Ad2jHVhTaCVLheCllByllOGBGOwn6Eu3Yex6tRKTUHlEUpSpRCbWiVB6FLZ1SykwgEi05LXEatIQkG2hepKXKTWo3toH249kYaK+XYBRcvr6XGwtjgBp6q3EBf/SESEp5Vkr5NFoi9R6wWi9fuVGprbr6ck6i/cA/QvEWU9BOKh65oYXOXmp16XfqxjhuuU+llOlSyteklIFo9af/KFLuctM2ies9d5T0mnYb8cWgJY4F/IF8IE5vIX9TStkM6AT0A/52G9tYkoIrAgOAk3rCCdq+XnLDvna6oXW9KVrZRUkSuZ58H0NLHm30mmIvIURDtCT5L5MyKWW2lPJTtBKfZnpsF26IzUVK2efWSypc3jIpZWe0/SvRjtebJruNRS1Hu0JSD60k67vSVllCDLuBPLQW7mcopdwDCuu1SzuWTpQ2n5TyqJQyXErpKaXshVY6tLfI+E+llEFSSh89dhvg+A3rFmjlMuMAL8BaP8HdB5TWOq4oVYpKqBWl8trALUo09Euyi4CPhBA+AEKIOkVqSF3QksMUod209697DUhKeRntsv27QrspriVaq/RSff3DhBDeemwFrZYl9QQQh36T0y0sQ/sB7wKsKjJ8IfC2nsAghPAWQgy4y02KA/yEELbw1/tUCNFPaDdOCiBN3zZTkWUV2yap1cuW9ipa7lOa5cCrQoj6ev1qQR1uvhDiIb2u2FqPxUjp+7pwG0uxAq32dyzFT16+QWu57iWEsNa/866i+M2t4Wg32JZkC1oZBGhlA53QEuHBaCVPnwOTpJT7SppZCDFBX5+Dnog/h3ZcH0JLCtOEdnOhgx5fsBCi7S22s2C5jYUQD+tXQHLQ/p2UtO8S0FpmSz1WpZSH9Ok+B36RUpbWWh8HBAit9ruor9GS1XwpZaldOEqtq8PSjqXmt9jWlvr35iiEeB3tBOdLfZy9vs+EEMIf+AztHohrNyxmBHBIaiVeSYCDEKIZWm1/uXZ1qSiWohJqRam8PgOG6slbaSajlSDs1ss6NqO1SgP8B3BAayXcjVa6cD88jVbrGQN8D/xLr60F6I1281IG2g2KQ6SUOSUsYx5aq941IcTHJYwHLZnsCmyR2s2QReddh1Z2kY62be3vclu2oPWocVUIUbCOW+3TIP3vDLSu4P4r9X6G0ep9/6mXH7x+l/HcaDFaq+U2tC7NctC6MAOt15XVaMl0BFqCWlI/xyVtYzFS68lhF1rCu7LI8MtordbT0JLGy8BE9N8WPXnNlFr3eSUtdzPgIYQYKqW8rJcO1ZZSviClfAjorpcdlCYb+ACtXCEReBl4XGo15Sa0qwSh+r5JREtq3W6xvAJ2aF3CJerL9tG38cb4s9Ba0P/Uv9cOpSxvOdCd4icjNyo4KUwSQhwsMnwJ2k2ApbZO36NngVi0WupuQI8iZSX2aDFnoJ2g7AKKPfBHCOEFjC8YLrX7N15BO64Wcv14VJQqTcjS75dRFKWCE0IsQ6urvJ9daVUrQogooKuUMsrCoVQ5QojvgP+TUm64xTR1gE1oN9stQmvRrI2WiIVIKR8pj1grKqHdDBgPtJJSnv2r6RVFsYyyfvCDoihlSEr5jKVjUJTSSCkfv41prgjtKZuvAWsBP7SW4bVc7x2mOhsL7FPJtKJUbKqFWlGUak0IMQH48ha1rYpiEfrVE4HW+496AqSiVGAqoVYURVEURVGUe6BuSlQURVEURVGUe1Cpa6i9vLxkQECARdadmZmJk1NJ3ecqSsWmjl2lslLHrlJZqWO3ajhw4ECilNK7pHGVOqEOCAhg//79Fln31q1b6dq1q0XWrZSP+HStNzcfF3sLR3J/qWNXqazUsatUVurYrRqEEBdLG1epE2pFKUt/X6bdA7RydEcLR6IoiqIoSkWmEmpFKcXYrg0sHYKiKIqiKJWASqgVpRRdG/tYOgRFURRFUSqBKpdQG41GoqOjyckp6WnG94+bmxsRERFlug6lZPb29vj5+WEwGMp0PTEp2QD4ujuU6XoURVEURancqlxCHR0djYuLCwEBAQghymw96enpuLi4lNnylZJJKUlKSiI6Opr69euX6bpeXXkYUDXUiqIoiqLcWpVLqHNycso8mVYsRwiBp6cnCQkJZb6uvz8cVObrUBRFURSl8qtyCTWgkukqrry+385BXuWyHkVRFEVRKjf1pERFKcWlpCwuJWVZOgxFURRFUSo4lVCXgXnz5hEcHEzz5s35z3/+Uzj8yJEjdOzYkRYtWvDoo4+SlpZW4vwfffQRzZs3Jzg4mKeffrrwBsvJkyfTsmVL/va3vxVOu2TJEubNm1e2G1SKTp06FX6eOHEizZs3Z+LEiSxcuJCvv/661PliYmJ44oknADh8+DAbNmwo81jvxsTVR5i4+oilw1AURVEUpYJTCfV9dvz4cRYtWsTevXs5cuQIP/30E2fPngVgxIgRzJkzh2PHjjFw4EDmzp170/xXrlzh448/Zv/+/Rw/fhyTycSKFStITU1l586dHD16FJPJxLFjx8jOzubLL7/kpZdeKu/NBGDnzp2Fn//3v/9x8OBB5s6dy5gxY4ol/Tfy9fVl9erVQMVOqF/t0YhXezSydBiKoiiKolRwKqG+zyIiIujQoQOOjo7Y2NgQHh7O999/D8Dp06fp0qULAD169OC7774rcRn5+flkZ2eTn59PVlYWvr6+WFlZkZeXh5SS7OxsDAYDc+fOZdy4caV2HxcbG0uXLl0IDQ0lODiY7du3A+Ds7Mxrr71Gq1at6NatW+ENfpGRkfTu3ZvWrVvz4IMPcurUKQDi4uIYOHAgISEhhISEFCbSzs7OAPTv35/MzEzat2/PypUrmTlzJu+//z4A586do3v37oSEhNCqVSsiIyOJiooiODiYvLw8ZsyYwcqVKwkNDWXlypUEBQUVxmM2m2nYsCGJiYn3/L3cjQ6BnnQI9LTIuhVFURRFqTzK7KZEIcRioB8QL6UM1ofVAFYCAUAU8JSU8po+biowHDAB46SUv9xrDG/+eIKTMSWXVdytZr6u/OvR5qWODw4OZvr06SQlJeHg4MCGDRto06ZN4bh169YxYMAAVq1axeXLl2+av06dOrz++uv4+/vj4OBAz5496dmzJwCPP/44YWFhdOvWDTc3N/bt28eMGTNKjWXZsmX06tWL6dOnYzKZyMrS6oEzMzNp1aoVH3zwAbNmzeLNN99k/vz5jBo1ioULFxIUFMSePXt46aWX2LJlC+PGjSs8MTCZTGRkZBRbz7p163B2dubwYa2buZkzZxaOGzp0KFOmTGHgwIHk5ORgNpuJj48HwNbWllmzZrF//37mz58PwKlTp1i6dCkTJkxg8+bNhISE4OVlmZsDIxO07Wzg7WyR9SuKoiiKUjmUZQv1l0DvG4ZNAX6TUgYBv+l/I4RoBgwBmuvz/FcIYV2GsZWZpk2bMnnyZHr06EHv3r0JCQnBxkY7b1m8eDGffvoprVu3Jj09HVtb25vmv3btGj/88AMXLlwgJiaGzMxMvvnmGwAmTZrE4cOH+eCDD3jjjTeYNWsWn3/+OU899RSzZ8++aVlt27bliy++YObMmRw7dqyw32wrKysGDx4MwLBhw9ixYwcZGRns3LmTJ598ktDQUEaPHk1sbCwAW7ZsYezYsQBYW1vj5uZ2W/siPT2dK1euMHDgQEB7IIujo+Mt53nxxRcL668XL17MCy+8cFvrKgvT1hxj2ppjFlu/oiiKoiiVQ5m1UEsptwkhAm4YPADoqn/+CtgKTNaHr5BS5gIXhBDngHbArnuJ4VYtyWVp+PDhDB8+HIBp06bh5+cHQJMmTdi0aRMAZ86cYf369TfNu3nzZurXr4+3tzcAgwYNYufOnQwbNqxwmkOHDgHQqFEjxo8fz7Zt2xgyZAhnz54lKOh638ldunRh27ZtrF+/nmeffZaJEyeWWNsshMBsNuPu7l7Yynw/SCnveJ66detSs2ZNtmzZwp49e1i6dOl9i+dOTerd2GLrVhRFURSl8ijvfqhrSiljAaSUsUIIH314HWB3kemi9WE3EUKMAkYB1KxZk61btxYb7+bmRnp6+n0O+2Ymk6nU9SQkJODt7c3ly5dZvXo1mzdvJj09vXC42WzmX//6F88///xNy/D09GTnzp3ExcXh4ODAxo0bCQsLKzbd1KlT+fjjj0lOTiYvL4/09HRMJhMJCQnUqlWrcLpLly7h6+vLkCFDSEpKYvfu3QwcOBCz2cySJUt44okn+OKLL2jXrh1CCPz9/fn6668ZOHAgUkqOHz9OixYt6NKlCx999BEvv/wyJpOJzMxMXF1dAYrFVfA5NzcXg8GAEILatWuzfPly+vXrR25ubmHJiNlsJj09HRsbG5KTk4stZ+jQoQwdOpQhQ4YUlqncKCcn56bvvqxsvVAuqyk3GRkZ5bbvFOV+UseuUlmpY7fqqygPdinpSR0lNm9KKT8DPgNo06aN7Nq1a7HxERER5fJI8Fs9erxPnz4kJSVhMBhYsGAB/v7+wPWSD9BanseOHYsQgpiYGEaMGMGGDRt4+OGHeeqppwgPD8fGxoawsDDGjRuHnZ0dAGvXrqVjx440aqT1PtG5c2c6depEy5Yti3VjB7B//37mzp2LwWDA2dmZr7/+GhcXF5ycnDh//jxdu3bFzc2NlStX4uLiwooVKxg7diwffPABRqORIUOG0KlTJ/773/8yatQoli5dirW1NQsWLKBjR+1x3EX3QcFnOzs77OzscHFxYdmyZYwePZp3330Xg8HAqlWrcHZ2xsrKChcXF/r06cO8efN48MEHmTp1KoMHD2bw4MG89NJLjB49utR9bG9vT1hY2N1+fbfl9FUtyW9cq2o9Yn7r1q3c+O9GUSoDdewqlVV1OHaNJjNZuSayjSZyjNp7wWftZSY7z0ROvonsPBO5+frfRabNNZqLzZ9jNJNjNJFr1KbXXiaWj+xAm4Aalt7kYsTdXJa/7YVrJR8/Fbkp8TTQVW+drg1slVI21m9IREr5rj7dL8BMKeUtSz7atGkj9+/fX2xYREQETZs2ve/bcqNbJdQVnbOz8003FlYk+/fv59VXXy3slaQk5fE9D/6fdvitHN2xTNdT3qrDf+xK1aSOXaWyqojHbm6+ibTsfNJzjGTk5pOZayIrL5+M3Hyy8kxkFhmWmZdPVq7p+ri8/OLjc03kmcx3HIOVAAeDNfaFLyscbK2xt7HGwdYaOxttmL3BGjsbK+xsrLEzWPFMO3/q1rj1PVllQQhxQErZpqRx5d1CvQ54Dpijv/9QZPgyIcSHgC8QBOwt59iUCmDOnDksWLDAorXTBab1KfsTM0VRFEW5U1JKso0mrmUZSc0ykpZjJD0nn7RsI+k5RtJytEQ5PSf/hnHa32k5+eTl314CbGdjhZOdDU521jjZ2uBoa42znQ01Xexx1Ic52dngZKslwY62NlpibLDGvkhyXDisSPJsa22FECUUKUgJeZmQmQCZiZAVB1lJ119WI4HyT6hvpSy7zVuOdgOilxAiGvgXWiL9rRBiOHAJeBJASnlCCPEtcBLIB16WUprKKrbqriK3Tk+ZMoUpU6ZYOgwAQuq6WzoERVEUpYozmSUpWXlacpydx7VMI9ey8kjJ0t6vZRn18cWH/VVC7GCwxsXeBlcHAy72Nrg72lK3hmPh3672BlztbXCxN+BsZ1OYNDva6smznQ2OBmtsrO9Th3CmfC0Zzkwo5ZUIGfHae2YC5GeXvBxrW2jUG9zr3p+47pOy7OXj6VJGdStl+reBt8sqHkW5UydiUgFo7nt73QQqiqIoCmjlFEkZeSRl5JGYkcuOaCMnt54r/DsxI7fwc3JmHuZSqm9trATujrZ4OBpwdzRQt4YjLf3c8HC0LTbc1d6Ai70BVwctQXaxt8FwvxLhvyIlZF+D9FjtlRYL6VchPUZ/14dlxoMs4STAygacvK+/vBqBs3fxYY6e11+2TlBSq7aFVZSbEhWlwpn140mg6tVQK4qiKHcnL99MQkYuV1NziEvLuf6uf05I15LltJz8m2c+fhonW2s8ne3wcrbFv4YjYf4eeDvbUsPJFg+n60myljBrLccllkSUp9x0SLkMqZch5ZL2nhoNaTHaK/0qmHJvns+hBrjUBtfaULO59tm55vUk2dkHnLzA3r1CJsh3SiXUilKKGY82s3QIiqIoSjnJN5mJTc3hSko20deyiUnJ5mpaDnGpWsIcl5ZDUmYeN/blYGttRU03O2q52tO0titezrZ4Odvh5WKHl7Mdns62RB4/RN/uXXC0rYBpV/Y1SL4AKRf1xDlaT54vQ+olyEktPr21LbjWATc/qNseXGqBq6/27lJbe3euBQb7+xKelJKknCQupF7gYtpFolKjeLbZs9R0qnlfln+/VMBvVlEqBlXqoSiKUnXk5puIScnhyrVsrqRkEX0tmyvXtOT5Sko2sanZN5VeeDgaqOlqTy03e1rUcSv8XMvVvvCzh6PhL1uR085bWS6ZllKrSU6+AMnntde1Ip+zrxWf3s4V3OpqCbN/e+2ze11w89fenXzA6v6Xk0gpicuK43zKeSJTI4lMieR86nkiUyJJy0u7Hp61HeF1w1VCXR3MmzePRYsWIaVk5MiRTJgwAYAjR44wZswYMjIyCAgIYOnSpYUPSCkqICAAFxcXrK2tsbGxoaBrwMmTJ/Pzzz8TGhpa+HjuJUuWkJyczPjx48tvA3WdOnVi586dAEycOJENGzbQp08fGjRogKOjY4lPZQSIiYlh3LhxrF69msOHDxMTE0OfPn3KM/TbcuRyCqBuTlQURakssvNMXEzOJCoxi4tJmUQlXf8cm5ZTrHXZSkAtV3v8PBxpV78Gfh4O1HF3oI6HA34ejtR2s8feYG25jblTeZmQeFZ/nYbEM3rSfAHyinRGIKy0ZLlGIDQfCB71tc8e9bTk2aHsf/MSsxM5nXyas9fOEpkaWZhEZxozC6dxt3OngXsDegf0JtA9kPqu9QlwC6CWUy2sRDnVh98BlVDfZ8ePH2fRokXs3bsXW1tbevfuTd++fQkKCmLEiBG8//77hIeHs3jxYubOnctbb71V4nJ+//13vLy8Cv9OTU1l586dHD16lKFDh3Ls2DEaNmzIl19+ycaNG8tr84opSKYB/ve//5GQkFD4AJpb8fX1ZfXq1QAcPnyY/fv3V8iE+p0NEYCqoVYURalIjCYzF5OyOBefwYXETKIStcT5YlIWV9Nyik1bw8mWep6OdAj0pG4NR+rWcKSOuwN+Hg7UcrMvvxv37qfMREjQE+bEM/rns1p5RgFhpSXKng2g3gN6wqwnzu7+YGNbLqHmm/O5mHaRU8mnOH3tNKeTtVdSTlLhNN4O3gS6BdK/QX8auDUg0D2QBu4NqGFfsR7c8ldUQn2fRURE0KFDBxwdtf4Rw8PD+f7775k0aRKnT5+mS5cuAPTo0YNevXqVmlDfyMrKiry8PK3vyexsDAYDc+fOZdy4cRgMhhLniY2NZfDgwaSlpZGfn8+CBQt48MEHcXZ2ZvTo0fz+++94eHiwYsUKvL29iYyM5OWXXyYhIQFHR0cWLVpEkyZNiIuLY8yYMZw/fx6ABQsW0KlTp8IHxPTv35/MzEzat2/P1KlTiYiIwNnZmddff51z584xZswYEhISsLa2ZtWqVVhbW9OvXz8OHjzIjBkzyM7OZseOHUydOpV//vOf7Ny5s/AR7Y0aNWL37t3FTi7Ky6wBweW+TkVRFEWTYzRxITGTs/EZnItL51xCBmfjMohKysRout7U7OVsS4CnEw809CLA05F6Xk7U93TC39MRN4eSfx8rhfxcSDgFcSe019Vj2ntW4vVpbBzAK0grzaFfKiUAACAASURBVPB6Vushw7uxljjb/HUD1/1kNBk5c+0MxxOPczL5JKeTT3Mu5Ry5+g2LBisDDd0b8qDfgzT2aEzjGo1p5NEIN7uqUV5ZtRPqn6doB+D9VKsFPDKn1NHBwcFMnz6dpKQkHBwc2LBhA23atCkct27dOgYMGMCqVau4fPlyicsQQtCzZ0+EEIwePZpRo0bh4uLC448/TlhYGN26dcPNzY19+/YxY8aMUmNZtmwZvXr1Yvr06ZhMJrKysgDIzMykVatWfPDBB8yaNYs333yT+fPnM2rUKBYuXEhQUBB79uzhpZdeYsuWLYwbN67wxMBkMt3Uj/W6detwdnbm8OHDAMycObNw3NChQ5kyZQoDBw4kJycHs9lMfHw8ALa2tsyaNYv9+/czf/58AE6dOsXSpUuZMGECmzdvJiQkxCLJNFS9R44riqJURCaz5GJSJhGx6UTEpnHqajrn4tO5lJxVWNNsJcC/hiMNfVzo1rQmQT7ONPRxJtDbCRf7Spw0F0i/quUrBUlz3Amt9bngkRw29uDTFBr3Bp9m4NVYS6Td6pZJPfNfMUszUalRHE86zrGEY5xIOsGp5FMYzUZAK9doXKMxQxoPoXENLXmu71Yfg1UV+K5KUbUTagto2rQpkydPpkePHjg7OxMSEoKNjbabFy9ezLhx45g1axb9+/fH1rbkSy5//vknvr6+xMfH06NHD5o0aUKXLl2YNGkSkyZNAmDEiBHMmjWLzz//nE2bNtGyZUv++c9/FltO27ZtefHFFzEajTz22GOEhoYCWmv34MGDARg2bBiDBg0iIyODnTt38uSTTxbOn5urnVVu2bKlsGbb2toaN7fbO5tMT0/nypUrDBw4EAB7+7++4/fFF19kwIABTJgwgcWLF/PCCy/c1rrKwoGLyQC0rle5LjspiqJUVBm5+ZyKTSMiNo2TegJ9+mo62UYtcbS2EgR6OdHc140BoXVo6ONMUE1nAjydKlc9861kJEDMoeKvjKvXx7vVhZrB0KSv1t1czWCtdMPKctufmpvK4fjDHIo/xLHEY5xMOkmGUWtcc7RxpJlnM4Y2HUqwVzDBXsH4Ovlavru/cla1E+pbtCSXpeHDhzN8+HAApk2bhp+fHwBNmjRh06ZNAJw5c4b169eXOL+vry8APj4+DBw4kL179xaWigAcOnQIgEaNGjF+/Hi2bdvGkCFDOHv2LEFBQYXTdenShW3btrF+/XqeffZZJk6cWOKNgkIIzGYz7u7uha3M94O8sW+h21C3bl1q1qzJli1b2LNnj0UfQf7vjacBVUOtKIpyNzJy8zkWncrR6BSORqdy7Eoql5KzCse7OxpoWsuVp9v506S2C81qu9LQx7nqJM4A2SkQcxD/i9/BikUQcxjSovWRQivRCOwKvmHaFfCazcvlpsBbkVISnR7NoYRDHIw7yOH4w0SmRgJgI2xoXKMxfQP7asmzZzD13epjbcFkv6Ko2gm1hcTHx+Pj48OlS5dYs2YNu3btKjbcbDYze/ZsxowZc9O8mZmZmM1mXFxcyMzMZNOmTTeVdbzxxht89tlnGI1GTCbtrN7KyqqwpKPAxYsXqVOnDiNHjiQzM5ODBw/yt7/9DbPZzOrVqxkyZAjLli2jc+fOuLq6Ur9+fVatWsWTTz6JlJKjR48SEhJCt27dWLBgARMmTMBkMpGZmVli7yQ3cnV1xc/Pj7Vr1/LYY4+Rm5tbGG8BFxcX0tPTiw0bMWIEw4YN49lnn8Xa2nL/SN8Z1MJi61YURalMcowmTsamcfSyljwfiU7hfGJmYa8adWs40KKOG4Pb1qVpbRea1nallqt91WrFlFLrVePyHu11aY9WA40kEKBGA/DvoCXPvmFQuyXYWb600CzNnL12lr1X93Io/hCH4g+RmK3VabsYXAj1CaVvYF/CfMII9grG3ub+9C9d1aiEugw8/vjjJCUlYTAY+PTTT/Hw8ABg+fLlfPrppwAMGjSosJwhJiaGESNGsGHDBuLi4gpLJPLz83nmmWfo3bt34bLXrl1L27ZtC1uxO3bsSIsWLWjZsiUhISHF4ti6dStz587FYDDg7OxcWLbh5OTEiRMnaN26NW5ubqxcuRKApUuXMnbsWGbPno3RaGTIkCGEhIQwb948Ro0axf/93/9hbW3NggUL6Njx9lptlyxZwujRo5kxYwYGg4FVq1ZhVaTe66GHHmLOnDmEhoYydepUBg8eTP/+/XnhhRcsWu4B0MDb2aLrVxRFqahiU7PZH3WN/VHJHLh0jVOx6eTrBc/eLnaE+GklGy393Gjp504Np/LpVaJcGXMg5qCeQO/V3rP03ivs3cCvHQQ/DnXbsiMyg849+lk2Xp2UkotpF9l7dS97Yvew7+o+ruVqfVHXca5Dh9odCPMJI8wnjAbuDSpkF3UVkbiby/IVRZs2bWRBH80FIiIiaNq0aZmvOz09HRcXy59Z3o2C3jkqqv379/Pqq6+yffv2Uqcpj+9593ntP8YOgZ5lup7ytnXrVrp27WrpMBTljqlj1zJMZsmpq2kcuHiN/VHXOHDxGldSsgFwMFgT5u9OmL87Leq4E1LXreq1PBfIz4Xo/RC1HaJ2aEl0wSO3PRtqTw2s205792pc7GZBSx+7VzOvsjt2N3tj97Ln6h7is7TOAWo61qR97fa0q9WOdrXaUdu5tsVirAyEEAeklG1KGqdaqJUKZc6cOSxYsMCitdMFPvr1DKBqqBVFqV7yTWaOXUll1/kkdkUmcehSChm5+QDUdLWjTb0aDO9cn7YBNWha2wWbytiX8+3Iz9NaoC9sh6htWgKdnwMIrd653Uio10lLoJ0s0xtVaYwmI4fiD7Hjyg62X9nOuZRzAHjYedCutpY8t6/dHn8X/6p58mMBKqGuhipy6/SUKVOYMmWKpcMAYO4TIX89kaIoSiVnNksirqaxKzKJnZFJ7L2QXJhAN6rpzGNhvrSpV4PW9Tzw83CougmYlJB0Ds79Buc2w8U/wajfm1SzBbR+Aeo/qCXRDh6WjbUEsRmxbL+ynR1XdrAndg9Z+VnYWNnQ2qc1A1oPoKNvR4I8glQJRxlRCbWilMLf09HSISiKopSJy8lZbD2TwJ9nE9l9IYmULK3/4EAvJ/qH+tKpgScdAj3xci7fh4OUu5xUuLBNT6J/u/60Qc+GEDYMAh6EgM7gWPG6TzVLMyeTTrLl0hZ+v/x7YSt0bafa9A3sS+c6nWlfuz1OBicLR1o9qIRaUUqx46x2l3PnoIp1KU9RFOVO5RhN7D6fxNbTCWw7k8D5xEwA6rg70KNpTTo19KRjoBe13Kp4Dw5Saj1vnN4AZzdrNxJKE9g6Q/1w6DwBGnYDjwBLR1oio8nIvrh9hUl0fFY81sKa1jVb83qb1+lcpzOBboFV9ypCBaYSakUpxSdbzgIqoVYUpXI6n5DB1tMJ/HEmgd3nk8jNN2NnY0XHBp78rWM9whv7EODpWPWTL1M+XNqlJdGnN8C1KG14rZbwwHgtgfZrBzYVsyeSLGMW269sZ8ulLWyP3k66MR0HGwce8H2Ah/0fpotflyrz+O7KTCXUilKKjwaHWjoERVGU22Y2Sw5dvsamk3H8eiKusBU60NuJoe3rEd7Ym/b1a1StB6eUJjddq4M+/TOc+QVyUsDaDgLDtSS6UW9w9bV0lKXKNeWyI3oHG6M28kf0H2TnZ+Nh50H3et152P9hOtTuoPqDrmBUQl0G5s2bx6JFi5BSMnLkSCZMmADAkSNHGDNmDBkZGQQEBLB06VJcXV35888/GTt2LHZ2dixfvpyGDRuSkpLC4MGD2bhxY7m3HsyYMYMuXbrQvXt3tm/fzpgxYzAYDKxfv57x48ezevXqUucdMWIE//jHP2jWrBnvvPMO06ZNK8fI7y9fdwdLh6AoinJLOUYTf55L5NeTcWyOiCMxIw8bK0HHBp48/0AADzX2oW6NanI/SE4anNkIJ77X6qFNueBQAxr3gcaPQIOHwa7iPl/AaDKyK3YXGy9sZMvlLWQaM/Gw86B/g/70CuhFK59W6omEFZjqh/ouldYP9fHjxxkyZAh79+7F1taW3r17s2DBAoKCgmjbti3vv/8+4eHhLF68mAsXLvDWW28xaNAg3nvvPaKioti4cSMffPABr732Gv379yc8PLzMt+VWxowZQ/v27e/qIStl2d91eXzPW09r/XR2bexTpuspb5buD1VR7pY6djWZuflsjojj52NX2XY2gaw8E852NnRt7E2PZjXp2tgHNweDpcMsH7npWgv0ie/h7K9aEu3iC80GQLP+WimHteXbDks7dqWUHIw/yI+RP/LrxV9Jy0vDxdaF7v7d6V2/N+1qtcPGyvLxKxrVD3U5ioiIoEOHDjg6ai0C4eHhfP/990yaNInTp0/TpUsXAHr06EGvXr146623MBgMZGdnk5WVhcFgIDIykitXrtwymZ4yZQrr1q3DxsaGnj178v777/P8889jb2/PiRMniIuL48MPP6Rfv36YTCamTJnC1q1byc3N5eWXX2b06NEA/Pvf/2bJkiVYWVnxyCOPMGfOHJ5//nn69etHSkoK3377Lb/88gubN2/m7bffpl+/fhw/fhyTycTkyZP55ZdfEEIwcuRI/v73v9O1a1fef/99Vq9eTXZ2NqGhoTRv3pzAwEC8vLwYP348ANOnT6dmzZqMGzeujL+Ru7dgayRQ9RJqRVEqn+w8E7+fjuenozFsORVPjtGMj4sdA8Pq0LN5LToE1sDOppq0XhqztVKOE2u0JDo/B1xqQ5sXofljWhJtVbG7hotOj+bHyB9ZF7mO6IxoHGwceNj/YR4JeIROvp0wWFeTE6IqpEon1O/tfY9Tyafu6zKb1GjC5HaTSx0fHBzM9OnTSUpKwsHBgQ0bNtCmTZvCcevWrWPAgAGsWrWKy5cvAzB16lRGjRqFg4MDS5Ys4fXXX+ett94qdR3Jycl8//33nDp1CiEEKSkpheOioqL4448/iIyM5KGHHuLcuXN8/fXXuLm5sW/fPnJzc3nggQfo2bMnp06dYu3atezZswdHR0eSk5OLrWfEiBHs2LGDfv368cQTTxAVFVU47rPPPuPChQscOnQIGxubm+adM2cO8+fP5/Dhw4VxDRo0iPHjx2M2m1mxYgV79+69vZ1uIZ88E2bpEBRFqcZyjCb+OJPAT0dj+S0ijqw8E17OtjzVpi59W9SmbUANrKyq+A2FBcxmuLQTjqyAkz9Abho414JWz0HzgdrDVSp4Ep1pzGRT1CbWRa5jf9x+BIJ2tdoxNnQs3f2742ioJqU5VVSVTqgtoWnTpkyePJkePXrg7OxMSEgINjbabl68eDHjxo1j1qxZ9O/fH1tb7Y7i0NBQdu/eDcC2bdvw9fVFSsngwYMxGAx88MEH1KxZs3Adrq6u2NvbM2LECPr27Uu/fv0Kxz311FNYWVkRFBREYGAgp06dYtOmTRw9erSw9jk1NZWzZ8+yefNmXnjhhcLW9Bo1br+fzc2bNzNmzJjCbfureQMCAvD09OTQoUPExcURFhaGp2fFfqS3j4u64UNRlPIlpWRf1DXWHIxm/dFY0nPz8XA0MCC0Do+2rE37QE+sq0sSDZB4Vkuij36r9RFt66yVc7QcrPUPXcFriqWUHEk4wpLEJUz6dhLZ+dnUc63H38P+zqOBj6pHfVchVTqhvlVLclkaPnw4w4cPB2DatGn4+fkB0KRJEzZt2gTAmTNnWL9+fbH5pJTMnj2blStX8sorr/Dmm28SFRXFxx9/zNtvv104nY2NDXv37uW3335jxYoVzJ8/ny1btgDcdAOjEAIpJZ988gm9evUqNu5ebniUUt7xvCNGjODLL7/k6tWrvPjii3e13vK0+WQcAN2b1fyLKRVFUe7NxaRM1hy8wppD0VxOzsbR1prewbUYEFqHTg08MVTVx3uXJCcVjq2Cw8vgygEQVhD4EHSbAU36gG3Ff1BJel46P53/iVVnVnH22lnshB39GvbjsYaPEeIdUvW7KqyGqnRCbSnx8fH4+Phw6dIl1qxZw65du4oNN5vNzJ49mzFjxhSb76uvvqJv3754eHiQlZWFlZUVVlZWZGVlFZsuIyODrKws+vTpQ4cOHWjYsGHhuFWrVvHcc89x4cIFzp8/T+PGjenVqxcLFizg4YcfxmAwcObMGerUqUPPnj2ZNWsWzzzzTGHJx+22Uvfs2ZOFCxfStWvXwpKPG+c1GAwYjUYMBq0WbODAgcyYMQOj0ciyZcvueL+Wt0XbzwMqoVYUpWyk5RjZcDSW7w5Gsy/qGkLAAw28eLV7I3oH18LRthr9REsJ0fvhwJdabbQxC3yaQ8/Z0OJJcKll6Qj/kpSSE0knWHVmFT9f+Jns/Gya1mjKjI4zcIl2oXen3pYOUSlD1ehfa/l5/PHHSUpKwmAw8Omnn+Lh4QHA8uXL+fTTTwEYNGhQsZ4zsrKy+OqrrwpbsP/xj3/w+OOPY2try/Lly4stPz09nQEDBpCTk4OUko8++qhwXOPGjQkPDycuLo6FCxcWloZERUXRqlUrpJR4e3uzdu1aevfuzeHDh2nTpg22trb06dOHd95557a2ccSIEZw5c4aWLVtiMBgYOXIkr7zySrFpRo0aRcuWLWnVqhVLly7F1taWhx56CHd3d6ytK/ZlOoAFw1pbOgRFUaoYKSWHL6ewfO8lfjwSS7bRRKC3ExN7NWZgWJ3q111n9jWtnOPAlxB/EgxOWgLd+jnwbQWVoCU315TLzxd+ZlnEMiKSI3CwceCR+o/wVKOnaO7VHICtMVstGqNS9lS3eXeptG7zLKmgd44nnnjC0qGUyGw206pVK1atWkVQUNA9Lau8vueqSHU9plRWlfnYTcsx8sOhKyzdc4lTV9NxtLWmf4gvg9vWJbSue/UrAYg+APsWad3d5eeAb5h2g2GLJ8CuYv22liYxO5GVp1fy7elvSc5JpqF7QwY3HkzfwL642Bbfhsp87CrXqW7zFIs7efIk/fr1Y+DAgfecTJeXjcdjAegdrG4aURTlzkkpORKdyrI9Fwtbo5vVdmX2Y8EMCPXFxb6adY2Wnwcn18Ke/8GV/WDrAqFDtdbo2iGWju62nUg6wdKTS/k56mfyzfmE+4UzrNkw2tdqX/1OjJRCKqGuQr788ktLh1CqZs2acf78eUuHcUe++DMKUAm1oih3Ji/fzPpjMXzxZxRHo1NxtLVmQKgvT7fzp6WfW/VLutLj4MAXsH8xZMSBZ0N4ZC6EPl1pWqPN0szWy1v56sRXHIw/iKONI081eopnmj5DPdd6lg5PqQCqZEJ9Nz1QKJVHeZUpLXquxKs6iqIoJUrMyGXp7kt8s+ciCem5BHo5MWtAcwaG1al+rdEAMYdg13+1sg6zEYJ6QrvR2iPAK3if0QWMZiMbzm9g8fHFnE89Tx3nOkxsM5GBQQNvKutQqrcql1Db29uTlJSEp6enSqqrICklSUlJ2NuXfR/RrtXxB1BRlDt2/EoqX/wZxY9HYsgzmQlv5M0LTwTQJci7+jx4pYCUELkF/pwHF/7QyjraDod2o8CzgaWju21Zxiy+P/c9X574kquZV2nk0Yj3HnyPngE91aPAlRJVuaPCz8+P6OhoEhISynQ9OTk55ZLUKTezt7cv7Nu7LP14JAaAR0N8y3xdiqJULlJK/jiTwMI/Itl9PhlHW2sGt63Lc50CaOjjbOnwyp8pX6uP/vM/cPWY9hTDHrOg9Qtg72rp6G5bam4qy08tZ2nEUlJyU2jl04o3OrzBg3UeVI10yi1VuYTaYDBQv379Ml/P1q1bCQtTj6auyr7ZfRFQCbWiKNflm8ysPxbLwj/OExGbRm03e6b1acLgtv64OVTDq1p5WXDoG9j1CaRcAq9G0H8+tHwKbOwsHd1tS81NZcnJJSyNWEqGMYNwv3CGtxhOmI/6nVduT5VLqBXlfvnyhXaWDkFRlAoix2hi1f7LfLb9PJeTs2ng7cTcJ1oyILQOtjaVox74vsrLhH2fw58fQ1Yi1G0Pvd+DRr0rTX00QFpeGt+c/IZvTn5DujGdHvV6MLrlaBrXaGzp0JRKRiXUilIKB9uK//AZRVHKVnqOka93XWTxjgskZeYR5u/OG32b0b1pzepXHw2Qm6H1H73zE8hK0m4w7DIJ6nW0dGR3JCMvg28ivuHrk1+TnpdON/9ujA0ZqxJp5a6phFpRSvH9oWgABoaVfb22oigVS3qOka92RrFo+wVSs410bezN2PAGtKtfo3rW0uamw97PYOd8yE6Ght0hfDLUrVxX8rLzs1kasZQvjn9BWl4aD9V9iJdCX6JJjSaWDk2p5FRCrSilWLH3MqASakWpTjJy8/VE+jwpWUa6N/VhfLdGtPBzs3RolpGXqT2IZefH2mPCG/aArlPAr3J1K5pvzmftubUsOLyA+Ox4wv3CeSn0JZp5NrN0aEoVoRJqRSnFNyPaWzoERVHKyY2JdLcmPozvHkRLP3dLh2YZJiMc/Br+eE97GEtQTwifAn6tLR3ZHZFSsuXSFv5z8D9EpUUR6h3K3PC5tKrZytKhKVWMSqgVpRQG68pzY42iKHcnN9/E0t2XmP/7OZIz83i4iQ/juwURUreaJtJmM5xYA1tmw7UL4N8Rnvoa/DtYOrI7diDuAB8e+JCjCUep71af/zz0Hx6u+3D1LNlRypxKqBWlFKv2ayUfT7apa+FIFEW530xmyQ+Hr/Dhr2eIvpbNAw09mdirCaHVNZGWEiJ/g81vwtWj4NMcnvlWa5muZAnopbRLvL//fX6//Ds+Dj7M7DiTAQ0HqAeyKGVKHV2KUorVB7SbElVCrShVh5SSrWcSeO/nU5y6mk5wHVfeHdSCB4O8LR2a5cQehU3T4cI2cPeHgZ9BiyfAqnL1dJSRl8Fnxz5jycklGKwMjAsbx7Bmw3CwcbB0aEo1oBJqRSnFytGVqxsoRVFu7fDlFN7dEMGeC8nU83Tk46fD6NeidvXs/g4gPQ5+nw0Hl4CDBzzyb+3Jhja2lo7sjpilmR/O/cC8g/NIykliQIMBjG81Hm/HanySpJQ7lVAriqIoVVpsajb/3nia7w9dwcvZjrcGNGdwW//q+UAWAGMO7P4vbP8A8nOg48vQ5XUtqa5kDsUfYs7eOZxMOkmIdwjzu80n2CvY0mEp1ZBKqBWlFMv3XgLg6Xb+Fo5EUZS7kZ1n4rNt51n4RyQmKXnloYaM7doAJ7tq+tMnJZz8AX59Q3tMeOM+0OMt8Gpo6cjuWGJ2Iu/vf5/159fj4+jDnAfn0Kd+H3XDoWIx1fR/FUX5az8djQFUQq0olY2Ukh+PxjJnQwQxqTn0bVGbKY80oW4NR0uHZjlxJ2HD63DxT+2Gw7/9AIFdLR3VHTOZTaw+s5p5B+eRbcpmZIuRjGgxAkdDNf5ulQpBJdSKUoqlIypfN1GKUt0djU7hzR9PcuDiNZr7uvLR4FDaB3paOizLyU2HrXNg9wKwd4V+H0Gr5yrdDYcAJ5NO8tautziedJz2tdozvcN06rvVt3RYigKohFpRFEWpAlKzjMzddIqley7h6WTLe4+34InWdbGurjccSgkn18LGqZAeC63+Bt1mglPlO7nIyMtg/uH5LD+1HHc7d1XeoVRIKqFWlFIs2RUFwLMdAywZhqIotyCl5PtDV3hnQwTJmXm80Kk+r/YIwsXeYOnQLCfxnFbecf53qNUCnloCddtaOqq7svniZt7Z8w6J2Yk81fgpxrUah6utq6XDUpSbqIRaUUqxOSIeUAm1olRUZ+PS+efa4+y5kEyYvztfvdiO5r5ulg7LcvJztZ47dnwENvZaN3hthoN15fupT8xO5J097/DrxV9pUqMJHz/8seq9Q6nQKt+/MkUpJ1+92M7SISiKUoKsvHw+/u0cn28/j5OdDe8OasHgNnWrb3/SAJf2wLq/Q+JpaPEk9HwbXGpaOqo7JqXkp/M/8d6+98g2ZjO+1Xiea/4cBqtqfMVBqRRUQq0oiqJUGtvOJDB1zTGupGTzRGs/pj7SBE9nO0uHZTm5GfDbLNj7Gbj5wdDVENTD0lHdlauZV3lz15vsuLKDUO9Q3nzgTQLdAi0dlqLcFpVQK0opFu+4AMCLndVd5IpiaalZRj4/lsuOjXsJ9Hbi29EdaVe/hqXDsqyzm+GnCZAaDe1GQbc3wM7F0lHdMSklq86s4sMDH2KWZqa0m8KQxkOwroQ9kSjVl0qoFaUUOyMTAZVQK4qlbTweyxs/nCApI5+XujZgXLcg7A3VONnKStZ67zi6Arwaw4u/gH97S0d1V+Iy45ixcwY7Y3bSvnZ7ZnaciZ+Ln6XDUpQ7phJqRSnF589VzrviFaWqiE/PYea6E2w4dpVmtV15OVjwfO8mlg7Lsk5v1Gqls5OhyyTtkeE2lbPk5ecLPzN792yMZiNvdHiDJxs9qbrCUyotlVAriqIoFYqUkrWHrzBz3UmyjSYm9mrMqC6B/Ll9m6VDs5ycNPhlGhxaAjWD4dk1Wpd4lVBqbipv736bn6N+pqV3S97p/A71XOtZOixFuScqoVaUUny2LRKAUV0aWDgSRak+kjPzmLbmGBtPXKV1PQ/ee7wlDX2cLR2WZV3YDmtfgrRo6PwP6Dql0rZK74zZyRt/vkFydjJ/D/s7Lwa/iI2VSkWUyk8dxYpSioMXUywdgqJUK79FxDH5u2OkZRuZ8kgTRj4YWH2fdAhgzNZ68Nj9X6gRqNVK162c3Xnm5Ofw0YGPWHZqGYFugXzy8Cc082xm6bAU5b5RCbWilGLhs60tHYKiVAsZufnM/ukkK/ZdpkktF5YMb0fT2tX8aXgxh2DNKEg8A21HQo83wdbJ0lHdlciUSCZum8jZa2cZ1nQY41uNx97G3tJhKcp9pRJqRVEUxWL2XkjmtVWHuXItm7FdGzChexB2NtW4Bw+zGXbN11qmnbxh2Bpo2M3SUd0VKSVrzq5hzt45OBocWdh9IQ/UecDSRan3zQAAIABJREFUYSlKmVAJtaKU4r9bzwHwUteGFo5EUaoeo8nMR7+eYcEfkdT1cOTb0R1pE1DN+5VOj4O1YyByCzTpB/0/AcfKuU/S89KZtWsWG6M20r52e97t/C7ejt6WDktRyoxKqBWlFCdj0iwdgqJUSZeTsxi34hCHLqUwuE1dZjzaDCe7av5z9P/s3Xd4VNXWx/HvTgiE3kvoxdA7oQpKl6tiAUF6k06oSlFEEaWoiAWkiUDovVoQCERESggl9N4JPbT0ZGa/fwze13uvg5FJsuck6/M89wkzKef33Ixn1jnZe60zW2BdP4h9BC9/CTW6g0VbyB2+fZgRO0ZwI/IGg6sPpkfFHngoD9OxhEhWRs5gSqmhQE9AA0eA7kAmYDlQHLgItNVa3zORTwiAaR2qm44gRKrz4+HrjFpzGDRMbV+NllUKmo5kVkKsY3nH7mmQrzx03Qj5yplO9VTs2k7AsQC+OfAN+TLlY36L+VTNV9V0LCFSRIoX1EqpQsAgoLzWOloptQJoB5QHArXWk5RSo4BRwMiUzieEECLpRcfZ+GjjMZbtu0K1ojn4pl01iuTKZDqWWXfPwcpucOOwY+Nh84/BK6PpVE/lYdxDRu8cTdCVIJoVa8bYemPJlj6NbywVaYqpv7GlAzIqpeJx3JkOA94FGj7+fAAQhBTUwqBvAs8AMKiJr+EkQljbiesPGbj0IOduR9C/YSmGNiuNl2caXwJwfD2sGwCe6aDdEij7kulET+1k+EmGbh/KjcgbjKo1ig5lO8jEQ5HmKK11yh9UqcHAeCAa2Ky17qiUuq+1zvGnr7mntc75F9/bG+gNkD9//hrLli1Lqdj/ISIigixZ0viwgVRuVmgMAH2qpK72TvLaFSlFa832KwksORlHZi9Fn8oZKJ/76Tt4pIbXrrInUPJ8AEWubuBhVl+OVRhJrLd1N+vtidjDivAVZPbITPc83SnpXdJ0JLeUGl67Aho1arRfa+33V59L8YJaKZUTWA28CdwHVgKrgGmJKaj/zM/PT4eEhCRnXKeCgoJo2LChkWML4Qp57YqUEBWXwLtrjrD+UBjPl87LF22rkCeLa9P9LP/afXANVnWHK3uhVh9o/gmkS2861VOJtcUyce9EVp9ZTe0Ctfn0uU/JnTG36Vhuy/KvXQGAUsppQW1iyUdT4ILW+jaAUmoNUA+4qZTy0VpfV0r5ALcMZBNCCOGis7ci6LdoP2dvR/BO89L0b/gMHml54iE4WuGt7unYhPjGXKjY2nSip3b10VWGBQ3jRPgJelXqxYCqA/D0SMO9w4XATEF9GaijlMqEY8lHEyAEiAS6ApMef1xvIJsQ/zZl8ykAhjUvYziJENbx4+HrjFgVSgYvTxb0qEUDX+suZ0gSdjvs+ByCJkLesvDmQshj3X0Zu8N2886v76C1ZmrjqTQs0tB0JCHcQooX1FrrvUqpVcABIAE4CMwGsgArlFJv4Si626R0NiH+LOxBjOkIQlhGvM3OxJ9OMvf3C1QrmoNvO1SnYA5rdqxIMjEPYHUvOPMLVG4HL0+x7PhwrTWLTixicshkSmYvyTeNvqFItiKmYwnhNox0+dBafwh8+F9Px+K4Wy2EW5jcporpCEJYwo0HMfgvOUDIpXt0q1ec914sR/p0abyLx+3TsKw93LsIL30Bfm9ZdlBLnC2OcbvHsf7cehoXacyEBhPI7GXNCwMhkksaH00lhBDCFSEXw+m76ABRcQl8074ar6T1QS0Ap3523Jn28nYMailWz3Sip3Y76jZDgoZw+PZh+lbpS78q/WTqoRB/QQpqIZz4dNNJAEa2KGs4iRDuaWnwZT5Yf5SCOTKypFdtSufPajqSWXY7/DYZto8Hn6rQbjFkL2w61VM7cvsIQ7YP4VH8I6Y0nEKzYs1MRxLCbUlBLYQT96PiTEcQwi3F2+x8/MNxFuy+RAPfPExtX40cmazZ/i3JxD6Cdf3gxEbHeumWX1l26iHAxnMbGbtrLHkz5WVh04WUySWbs4V4EimohXBiYqvKpiMI4XbuRsTSf/EB9l4Ip/dzJRnxQhnSpfWph+EXYGl7uHMKXpgAdfpbdr20XduZdnAa3x35jpoFavLF81+Q0/uJIyGEEEhBLYQQIpGOhT2g94L93I6I5cs3q/B6NesuZ0gyl3bDsg6g7dBpDZRqZDrRU4u1xfL+zvfZdHETrXxb8X6d9/Hy8DIdSwhLkIJaCCfG/3gcgNEvlTecRAjzfj5ynaErDpEzU3pW9a1L5cI5/v6bUrvQ5bDBH7IXgY4rIXcp04meWnhMOIO2DSL0dihDawyle4XuKIveZRfCBCmohXAiJt5uOoIQxmmtmfHrOT7bdIrqRXMws3MN8mX1Nh3LLLsdgiY4BrYUbwBtF0CmXKZTPbXz98/TP7A/d6LvyOZDIZ6SFNRCOPHxaxVNRxDCqLgEO6PXHmHl/qu0rFKQz9+ojLdXGh8xHR/t2Hx4bC1U6wQvfQnprLshc8/1PQzbPoz0numZ98I8KuWtZDqSEJYkBbUQQoj/cT8qjr6L9rPnfDiDmvgytKmvLAF4dNMxrOXaAWj6ETw72LKbDwHWnlnLuN3jKJ69ON82+ZaCWaSHuBBPSwpqIZz4aOMxAD5sWcFwEiFS1sU7kfSYv4+r96KZ0rYKrarL5kNuHoclbSHyDry5EMq1NJ3oqWmtmRE6gxmhM6hXsB6Tn59M1vRpvIe4EC6SgloIIcS/BV8Ip/fCEBSwqGdtapWw7trgJHPhN1jW0dFXusfPULCa6URPLcGewCd7PmH1mdW8WupVPqz3oXTyECIJSEEthBNyZ1qkNesOXmP4qlCK5MzE3G41KZ4ns+lI5h1dDWv7Qs4S0Gk15ChiOtFTi06IZsSOEQRdCaJXpV4MrDZQlvEIkUSkoBZCiDROa82sHeeZ9PNJapfIxazONWTyIcCuabB5NBSt5xgjbuFOHvdj7uO/zZ/Dtw/zXu33aF+2velIQqQqUlAL4cSYdUcB6fYhUje7XTPuh+PM33WRlyr7MKVtFTKkS+OdPOx2RyG9ZzqUfxVenw1e1m0VGBYRRt+tfbn26BpfNPxC2uIJkQykoBbCCW+vND5OWaR6MfE23l4Ryo9HrtPj2RK8/1I5PDzS+BKA+BhY19fRFq92X8cocQ/rXmCcCj9Fv639iEmIYVazWfgV8DMdSYhUSQpqIZyQCYkiNXsQHU/vBSHsvRDOey+WpVeDkrKeNvqeY/Phpd+h2cdQb6Cl2+KF3Ahh4LaBZPLKRMC/AvDN6Ws6khCplhTUQgiRxtx4EEO3ecGcvRXBV29W5bVqhUxHMu/hdVjUCu6cgVZzoHIb04lcsuPqDoYFDaNgloLMajoLnyw+piMJkapJQS2EE++uOQzAxFaVDScRIumcvfWILt8H8yA6nnnda9LAN6/pSOaFn4cFrzl6THdcCaUamU7kkk0XNvHub+/im9OXmc1mksvbupsphbAKKaiFcEK6HIjUZv+lcHrMD8HL04PlfepSsVB205HMu3kMFr4OtjjouhEK1zCdyCWrTq9i3O5xVMtXjWlNpsnAFiFSiBTUQjgxskVZ0xGESDK/nr5Nn4Uh+GTPSED3WhTNncl0JPOu7IPFbzgGtnT/GfKVM53IJfOPzueL/V9Qv1B9pjScQsZ0GU1HEiLNkIJaCCFSuZ+OXGfwsoP45svKgrdqkSdLBtORzDu3zbEBMUt+6LIechYzneipaa2ZenAq3x35jhbFWzCh/gS8PGX6oRApSQpqIZx4Z2UoAJPbVDGcRIint2LfFUatOUz1ojn5vltNsmeUQotj62B1T8hbBjqtgaz5TSd6anZtZ1LwJJaeXEpr39aMqTMGTwu3+RPCqqSgFsKJgtmtO8hBCIA5v53nkx9P0MA3D7M61yBTejnlc2ABbBwMhWtChxWQMYfpRE/NZrfxwa4P2HBuA90qdGNYjWHS+lAIQ+TsKoQTw5qXMR1BiKeitebLrWf4JvAML1YqwJdvVpXphwC7p8Mv78IzTaHtAkif2XSip5ZgT+C9397j54s/41/Vn96Ve0sxLYRBUlALIUQq8udR4m39CjOxVWU80/r0Q4DfpkDgR1DuFWj9PaSzbhefeHs8I3eMZMulLQyrMYzuFbubjiREmicFtRBODFl2EICv2lUznESIxLHZNSNXH2bV/qu8Vd8xSjzN37XUGoImwa+ToFIbeG0meFr3rS/OFsc7v77D9ivbGVFzBJ3LdzYdSQiBFNRCOFUybxbTEYRItASbnaErQtkYGsbQpqUZ1OQZKaa1dtyV3vklVO0Ir0wFC2/Yi7XFMixoGDuu7uC92u/Rvmx705GEEI9JQS2EE4Oa+JqOIESixNvsDFp6kJ+P3uDdf5Wlz/OlTEcyT2v4ZTTs+RZqdIeXpoCHh+lUTy0mIYYh24fwe9jvfFD3A9qUtvZodCFSGymohRDCwmITbPgvOciW4zcZ83J53qpfwnQk8+x2+OkdCPkeaveDFhPBwnfroxOiGbhtIMHXgxlXbxyv+75uOpIQ4r9IQS2EE/5LDgAwrUN1w0mE+Gsx8Tb6LdrP9lO3GfdqBbrULW46knl2m6Mt3sGF8OxgaPqRpYvpqPgoBgQO4MCtA4yvP56WpVqajiSE+AtSUAvhRPmC2UxHEMKp6DgbvReGsPPsHSa2qkT7WkVNRzLPboN1/eDwcnh+JDR81/LFdL+t/Qi9HcrE+hN5seSLpiMJIZyQgloIJ/o3fMZ0BCH+UlRcAm/ND2HPhbt81roybfyKmI5knt0G6/o7iunG78Nzw00ncskfd6ZDb4cy6blJtCjewnQkIcQTSEEthBAWEhGbQI95+wi5FM6XbavyWrVCpiOZZ7fB+gFweFmqKKajE6IZtG0QB24dYGL9iVJMC2EBUlAL4UTfhfsBmNm5huEkQjg8iomn27x9HLpyn6/bVaNllYKmI5lnt8OGgRC6FBqNtnwxHWuLZfC2wQTfCGZ8/fGyzEMIi5CCWggnqhfLYTqCEP8WEZtAt3n7CL1yn2ntq/GvSj6mI5lnt8PGQXBoMTw/Cp4fYTqRS+JscQzZPoTd13czrt442YAohIVIQS2EE72fk16+wj1ExTmWeRySYvr/2e3ww+NuHs+NgIajTCdySbwtnreD3mbntZ18WPdDaY0nhMVYt8u9EEKkAdFxNt6aH0LIpXC+erOqFNPgKKZ/HAYHFkCDt6HRe5bu5hFvj2f4juEEXQ1idO3RvFH6DdORhBD/kNyhFsKJngH7AJjTtabhJCKtiol3tMbbc+EuU9pWkTXT4JiA+NM7sH8e1B8KjcdYuphOsCcwascoAi8HMqrWKNqVbWc6khDiKUhBLYQT9UrlMR1BpGGxCTb6LNzPzrN3+Kx1ZV6vVth0JPO0hp+GOyYg1hsETT60dDFts9t4b+d7bL60mXf83qFjuY6mIwkhnpIU1EI40UNGOAtD4hLs9F90gF9P32ZSq0rSZxpAa0qdmwtXN0Bdf2g2ztLFtNaacXvG8fOFnxlcfTBdK3Q1HUkI4QJZQy2EEG4k3mZn4NIDBJ68xcevVaSdTEB0CJpIkasboFZvaP6J5Yvpz/Z9xpoza+hVqRc9K/U0HUkI4SIpqIVwouvcYLrODTYdQ6QhCTY7Q5Yd4pdjNxnbsjyd6xQzHck97PwKfv2U6wWaQotPLV1MA8wIncGiE4voWK4jA6sNNB1HCJEEZMmHEE40LZfPdASRhtjsmmErQvnxyHXef6kc3Z6VJUcABH8HWz+ECq04lacTPh7Wvg8UcCyAGaEzeO2Z1xhRcwTK4hcHQggHKaiFcKJz3eKmI4g0wm7XjFx9mA2hYYxsUZaeDUqajuQeDi1xdPQo8yK0mg2//W46kUtWnl7J5JDJNC/WnLF1x+KhrH1xIIT4f/JfsxBCGKS1ZtwPx1m1/yqDm/jSr6EMFALg2FpYPwBKNoI35oGnl+lELvnp/E98vPtjGhRqwKQGk/D08DQdSQiRhKSgFsKJjnP20HHOHtMxRCr3xebTzN91kZ71SzCkqa/pOO7h9C+wuicUqQ3tFoOXt+lELtl+eTvv7XyPGvlrMKXhFLwsfnEghPhfsuRDCCderixDNETymhF0jmnbz9KuZhFGv1RO1tMCnP8VlneGApWgw3JIn9l0IpfsDtvN27++Tfnc5ZnWZBre6ax9cSCE+GtSUAvhRHtpVyaS0cLdF/l000laVinI+NcrSTENcHkvLG0PuUtBpzXgnd10IpccunWIwdsHUzx7cWY0nUFmL2tfHAghnJMlH0IIkcLWHLjKmPXHaFouH1PaVsHTQ4pprh+GxW0gmw90XgeZcplO5JIz987QP7A/+TLlY3az2WTPYO2LAyHEk0lBLYQTb87azZuzdpuOIVKZTUdvMHzVYeqVys20DtXx8pTTMHfPwaLWkCGLo5jOmt90Ipdci7hG3y19yeiZkVnNZpEnYx7TkYQQyUyWfAjhxBs1CpuOIFKZHadvM2jpQSoXzs53Xfzw9pJODzy8DgtfA3sCdPsRclh7zPrd6Lv02dKHaFs0AS0CKJSlkOlIQogUIAW1EE608bP2G7twLyEXw+m9MIRS+bIwv1stMmeQ0y/R92BRK4gKh64bIG9p04lcEhkfSf/A/tyMvMns5rPxzSldW4RIK+SMLoQT8TY7gPxJXrjsxPWHdJ+/j4LZM7LwrVpkzyRt04iLhCVvwt2z0HElFKphOpFL4mxxDN42mFPhp/im8TdUy1fNdCQhRAqSgloIJzrN2QvA8j51DScRVnb5bhRd5gaTOX06FvasTZ4sGUxHMi8hDlZ0gav7oE0AlGxoOpFLbHYbo34bxd4be5lQfwLPFX7OdCQhRAqTgloIJ9rVkiUfwjW3HsXQee5e4m12lvSpS6EcGU1HMs9uh3X94OxWaPkNlH/FdCKXaK0Zv3c8Wy5tYbjfcFqWamk6khDCACmohXDi9WqyKVE8vYcx8XSbu49bD2NZ3Ks2vvmzmo5kntawaSQcXQVNx0KNrqYTuWx66HRWnl7JWxXfokuFLqbjCCEMkcWhQjgRHWcjOs5mOoawoJh4G70CQjhz6xEzO9egetGcpiO5h18/heDZUG8gPDvEdBqXLT6xmJmhM3n9mdcZXH2w6ThCCIPkDrUQTnSbFwzIGmrxzyTY7AxcepDgi+F83a4az5fOazqSe9g7G4ImQtVO0OxjsPhkyJ/O/8Sk4Ek0KtKID+p+IJMuhUjjpKAWwolOdYqZjiAsRmvNe2uPsOX4TT56pQKvVCloOpJ7OLwSfh4OZV6Cll9bvpjedW0Xo38fTY38Nfjsuc9I5yFvpUKkdXIWEMKJllIMiX9o0qaTrAi5yqAmvnStV9x0HPdwbhus6wvF6sMbc8HT2m87x+8eZ0jQEEpmL8k3jb/BO5236UhCCDcga6iFcOJhTDwPY+JNxxAWMXvHOWb9ep7OdYoxtKkM9AAg7BAs7wx5ykD7JeBl7eLz6qOr9N/anxwZcjCj6Qyypc9mOpIQwk0YKaiVUjmUUquUUieVUieUUnWVUrmUUluUUmcef5RdPMKoXgEh9AoIMR1DWMDKkCtM+OkkL1f2YewrFWQ9LcC9i7C4DWTMCZ1WgXd204lccj/mPv229iPOHseMpjPIlymf6UhCCDdi6g7118AmrXVZoApwAhgFBGqtfYHAx4+FMKb7s8Xp/mxx0zGEm9t+8haj1hyhgW8eprStiqeHFNNE3oWFrcAWB51WQzZrL5+KSYjBf5s/YRFhTG08lVI5SpmOJIRwMym+mE0plQ14DugGoLWOA+KUUq8CDR9/WQAQBIxM6XxC/KFFRR/TEYSbO3TlPv0XH6CcT1ZmdKpB+nSyis4xUrwtPLwGXdZD3jKmE7nEZrcxcsdIDt8+zOTnJ1Mjv7VHpAshkoeJs39J4DYwTyl1UCk1RymVGcivtb4O8Pij/D1NGBUeGUd4ZJzpGMJNXbwTSY/5+8iTNT1zu9UkSwZrb7ZLErYEWNUDwg5A6++haB3TiVyitWZi8ES2XdnGyFojaV68uelIQgg3pbTWKXtApfyAPcCzWuu9SqmvgYfAQK11jj993T2t9f+so1ZK9QZ6A+TPn7/GsmXLUij5f4qIiCBLlixGji1SxsS90QC8Wzt1jYuW167rHsRqPtkTTUyCZnSdjBTILHem0Zoyp6bhc2Mrp337ElboX0l+iJR+7W5+sJmN9zfSJFsTXsv5WoodV6Q+ct5NHRo1arRfa+33V58zcUvlKnBVa7338eNVONZL31RK+WitryulfIBbf/XNWuvZwGwAPz8/3bBhwxSI/L+CgoIwdWyRMhLy3QSgYfn8hpMkLXntuiYyNoF2s/cQkRDLkl61qSZTEB22T4AbW+G54ZRu/D6lk+EQKfna3XhuIxsvbeRfxf/FpOcm4aHkokk8PTnvpn4pXlBrrW8opa4opcporU8BTYDjj//XFZj0+OP6lM4mxJ81TWWFtHBdvM1Ov8UHOH79Id91qSHF9B9C5jnGilftBI1Gm07jst1hu/ng9w+oVaAWn9T/RIppIcTfSnRB/Xidc4zW2pYExx0ILFZKpQfOA91xrOdeoZR6C7gMtEmC4wjx1G49igEgX1Zr984VSUNrzcjVh9lx+jaTWlWicVm54ALg5E/w4zB4phm0/MryUxBPhp9kaNBQSuQowZeNviS9Z3rTkYQQFuC0oFZKeQDtgI5ATSAWyKCUug38BMzWWp95moNqrQ8Bf7UGpcnT/DwhksPAJQcBWN6nruEkwh1M3nyKNQeuMbRpadrVKmo6jnu4EuzYhOhTFdoGgKeX6UQuCYsIo//W/mTxysL0JtNlcIsQItGedId6O7AVeBc4qrW2AyilcgGNgElKqbVa60XJH1OIlNevofSaFQ4Ld1/k2+3naF+rKIOaPGM6jnu4fdrRHi+bD3RcCekzm07kkgexD+i3tR8xCTEs+NcCCmQuYDqSEMJCnlRQN9Va/8/cZa11OLAaWK2UsvbtCCGeoGEZ6dwoYNPRG3yw4RhNy+Xn41dlCiIAj27CotbgkQ46rYHMeUwnckmsLZZB2wZx5dEVZjWbxTM55aJJCPHPOC2o/7uYVkp5A52AjMASrfXdvyq4hUgtwu472uYVzJG62uaJxNt3MZxByw5StUgOpravRjpP2Zz278EtUXeg24+Qq4TpRC6xaztjdo7hwK0DfP7c59QsUNN0JCGEBf2Td4evAU8gBliXPHGEcB9Dlx9i6PJDpmMIQ87eekTPgBAK58jI911rkjG9p+lI5tltsOotuHEY3pgLhaqbTuSyaQen8fPFnxlSfQgtSrQwHUcIYVFP2pS4BBijtT73+KlcwOLH/x6c3MGEMG1gY1/TEYQhtx7F0HXuPrw8PQjoUYtcmaXTA1rDplFw+md4cTKUSfrBLSlt9enVfHfkO1r7tqZHxR6m4wghLOxJa6jfBz5RSoUBHwOTgQ2ANzA2+aMJYVZ9X2uvCxVPJyougZ4BIYRHxrG8Tx2K5MpkOpJ72DMDgmdDXX+o1ct0GpfturaLj/d8zLMFn2V0ndGyNl4I4ZInraE+D3RQStUHlgM/As2SqA+1EG7v8t0oAIrmloIqrbDZNYOWHuTotQfM7uxH5cI5TEdyDyc2wi/vQblXoNnHptO47PS90wz7dRglc5Rk8vOT8fKQ/fVCCNc4XUOtlMqplBoAlAfaAg+AX5RSL6dUOCFMGr4qlOGrQk3HEClEa824jcfYeuIWY1+pIJMy/3A1BFb3gsJ+0Go2eFh7Y+atqFsMCBxA5nSZmd5kOlnSZzEdSQiRCjzpzLgOxzAXb2Ch1noB0BKooZTakBLhhDBpaLPSDG1W2nQMkUK+33mBgN2X6Fm/BF3qFjcdxz2EX4Alb0LW/NBuKXhZu+NNVHwU/oH+PIh9wLQm06TXtBAiyTxpDXVuYAmONnldALTW0cBHSimfFMgmhFF1SuY2HUGkkE1HrzP+pxP8q2IB3nuxnOk47iEqHBa3AXsCdFwFWfKaTuSSBHsCw3cM59S9U0xtPJVyueX3LIRIOk8qqD8EtgA2YNSfP6G1vp6coYRwB+duRwBQKq/8STg1O3D5HoOXHaJqkRx8+WZVPDxkcxoJsbC8M9y/BJ3XQR5rd7zRWjMpeBI7ru5gTJ0xPFf4OdORhBCpzJM2Ja7GMRFRiDTpvTVHAFjep67hJCK5XLobSa+AEPJn82ZOFz+8vaTXNFrDen+4tBNazYHiz5pO5LIFxxew/NRyulXoRtsybU3HEUKkQk/qQz0b+EZrffQvPpcZeBOI1Vov/p9vFiIVGNGijOkIIhndi4yj+7x92LRmfvea5M6SwXQk97B9AhxZAY3HQOU2ptO4bOulrXwR8gXNijVjaI2hpuMIIVKpJy35mA58oJSqBBwFbuPYoOgLZAPm8v+DXoRIdWoUy2U6gkgmMfE2ei8M4eq9aBb3qk1JWdbjcHAR7PgMqnWGBm+bTuOyw7cPM+q3UVTKW4kJ9SfgoazdoUQI4b6etOTjENBWKZUF8AN8gGjghNb6VArlE8KYUzceAVCmQFbDSURSsts1w1cdZt/Fe0xtX42axeXCCYBz22HjYCjZCF7+Eiw+6OTKoysM3DaQvBnz8k2jb/BO5206khAiFXvSHWoAtNYRQFDyRxHCvXyw3rHaSdZQpy6fbz7FxtAwRrYoS8sqBU3HcQ83j8OKLpCnDLQNAE9rDzp5EPuAAYEDSLAnML3pdHJnlI49Qojk9bcFtRBplbRPS32W7L3MjKBztK9VlL7PlzQdxz08vO5oj+eVCTquAO/sphO5JM4Wx9CgoVx9dJXZzWZTInsJ05GEEGmAFNRCOFGliIydTk22n7rFmPVHeb50Xj5+tQLK4ksakkRsBCx9E6LvQY+fIXth04lcorVm7K6x7Luxj4kNJuJXwM90JCFEGvHEHRpKKU+l1OcpFUYId3Is7AHHwh6YjiGSwLFxrRJVAAAgAElEQVSwB/gvPkCZ/Fn5tmN10nnK5jRsCbCqB9w4Am3mg08V04lcNiN0BhvPb8S/qj8vl3zZdBwhRBryxHcVrbUNx6hxuZUj0pxxG48zbuNx0zGEi8LuR9Nj/j6yZfRibreaZMkgf5hDa9g0Es78Ai9OhtLNTSdy2fqz65kROoNXS71K78q9TccRQqQxiXlnOQisV0qtBCL/eFJrvSbZUgnhBj5oWd50BOGiRzHx9Ji/j8hYGyv71qVAdun0AMDuabBvDtQbBDXfMp3GZXuv72XsrrHU9qnNh3U/lOU8QogUl5iCOhdwF2j8p+c0IAW1SNUqFLT25qy0Lt5mp//iA5y9FcG87jUp55PNdCT3cGwdbH4fyr8KTT8yncZl5+6fY+j2oRTLVowpDafgZfEOJUIIa0pM27zuKRFECHcTeuU+IJsTrUhrzftrj/LbmTt81royDXzzmo7kHq4Ew9o+ULgWvD4LPKy9lvxO9B36b+1PhnQZmN50OtnSy0WTEMKMvz2bKqUKK6XWKqVuKaVuKqVWK6WsvRVciESY8NMJJvx0wnQM8RS+3X6W5SFXGNj4GdrWLGI6jnsIPw9L20FWH2i/FLwymk7kkuiEaAZtG8S92HtMazyNglmkp7gQwpzELPmYBywB2jx+3Onxc82SK5QQ7mDcqxVNRxBPYd3Ba0zefJrXqhZkWLPSpuO4h6hwR69pbYeOqyBzHtOJXGLXdt777T2O3jnKV42+okKeCqYjCSHSuMQU1Hm11vP+9Hi+UmpIcgUSwl3IyHHr2XP+LiNWHaZ2iVx8+kZl2ZwGEB8DyzrA/cvQZQPkecZ0Ipd9uf9Ltl7eyoiaI2hctPHff4MQQiSzxCygu6OU6vS4J7WnUqoTjk2KQqRq+y+Fs/9SuOkYIpHO3npE7wUhFMmVkdmd/ciQztN0JPPsdljfHy7vhtdmQLG6phO5bMWpFcw/Np92ZdrRqVwn03GEEAJIXEHdA2gL3ACuA288fk6IVO2zTaf4bNMp0zFEItx+FEu3eftIn86D+d1rkT2TdHoAYPsncHQ1NPkQKr1hOo3LjkcfZ8LeCTQo1ICRtUbKXyCEEG4jMV0+LgOvpEAWIdzKhFaVTEcQiRAVl0DPgH3ciYhlee+6FMmVyXQk97A/AH77Aqp3hfpDTadx2anwU8y7PY9ncjzD589/TjoPGdAjhHAfTs9ISqkRWuvPlFJTcfSd/g9a60HJmkwIw0rlzWI6gvgbNrtm8LJDHL72gFmdakiLwz+cDYQfhkKpJvDSF2DxO7m3om4xIHAAGTwyMK3JNDJ7ZTYdSQgh/sOTLvH/6BcWkhJBhHA3e847tgrUKZnbcBLhzMc/HGfL8Zt82LI8zSsUMB3HPdw4Ciu6Qr5y0GY+WHzQSVR8FP6B/jyMe8igfIMokFl+z0II9+O0oNZab1RKeQIVtdbDUzCTEG7hyy2nAVjex/obuVKjuTsvMH/XRXo8W4Luz5YwHcc9PAyDJW0hQxbosAK8rT3oxGa3MXLHSE7dO8XUxlOxn7WbjiSEEH/piYvQtNY2pVSNlAojhDv5/I0qpiMIJzYdvcHHPx7nhQr5Gf1SOdNx3EPsI0cxHfMAuv8M2QuZTuSyySGTCboaxHu13+O5ws8RdDbIdCQhhPhLidnVcVAptQFYCUT+8aTWek2ypRLCDRTNLZvb3NHBy/cYsvwglQvn4Ks3q+HpYe31wUnClgAru8PN44470z6VTSdy2eITi1l0YhGdynWifdn2puMIIcQTJaagzoWj7/Sfu+drQApqkartPHMHgPq+1p4ql5pcvhtFz4AQ8mbNwPdd/ciYXnpNozX8PBzOboGXvwLfpqYTuezXK7/y2b7PaFikIe/4vWM6jhBC/K3EtM3rnhJBhHA3U7edAaSgdhf3o+LoNj+YBLtmfvda5MmSwXQk97DrGwiZ62iN52f90/WJuycYvmM4ZXOV5dMGn+LpIRdNQgj397cFtVKqNDADyK+1rqiUqgy8orX+JNnTCWHQl29WNR1BPBabYKP3wv1cDY9m4Vu1pKXhH46thS0fQIVW0PgD02lcdiPyBv6B/mTPkJ1pjaeRyUuWXQkhrCExkxK/A94F4gG01oeBdskZSgh3UDBHRgrmyGg6Rppnt2uGrzxM8IVwPm9TmdrSxtDh8h5Y0weK1HGMFfdIzOncfUXGRzIgcACRCZF82+Rb8mbKazqSEEIkWmLWUGfSWgf/14jXhGTKI4TbCDp1C4CGZfIZTpK2fbHlFBtCwxj+QhlerWr9zhVJ4u45WNoesheG9kvBy9t0Ipck2BN459d3OHf/HN82+ZbSOUubjiSEEP9IYgrqO0qpUjyelqiUegO4nqyphHADM4LOAVJQm7Qs+DLfbj9Hu5pF6N+wlOk47iHyDixq7Zh+2GkVZMplOpFLtNZMCp7Ezms7+aDuBzxb6FnTkYQQ4h9LTEE9AJgNlFVKXQMuAB2TNZUQbmBqh2qmI6Rp20/dYvS6ozxXOi8fv1YRZfHx2UkiPhqWtoNH16HrD5CrpOlELlt4fCHLTy2ne4XutCndxnQcIYR4KokpqLXWuqlSKjPgobV+pJSSsWQi1cuX1dp/Rreyo9ceMGDxAcrkz8r0jtXx8rT2+uAkYbfDmt5wNQTaLoAiNU0nclng5UAmh0ymWbFmDKkxxHQcIYR4aol5l1oNoLWO1Fo/evzcquSLJIR72Hr8JluP3zQdI825ei+K7vP3kSOjF/O61yRLhsRc96cBW8bAiQ3wwngo/4rpNC47eucoo3aMolKeSkyoPwEPJRdNQgjrcvpOpZQqC1QAsiulWv3pU9kAuXUnUr3vfjsPQNPy+Q0nSTseRMfTfd4+YuJsLOpXj/zZ5FQDwN7ZsHsa1OoDdfqbTuOysIgw/AP9yZ0xN183/hrvdPJ7FkJY25Nu/ZQBXgZyAC3/9PwjoFdyhhLCHczoVMN0hDQlLsFO34X7uXg3koDutShTIKvpSO7h5E+waSSUeQlaTHRsRrSwR3GPGBA4gDhbHN+/8D15MsrgJCGE9TktqLXW64H1Sqm6WuvdKZhJCLeQK3N60xHSDK01I1cfZvf5u0xpW4V6z0iRBcC1/bCqB/hUhdZzwOJTA+Pt8bwd9DYXH1xkRrMZlMohnVuEEKlDYhatva6UyqaU8lJKBSql7iilOiV7MiEM23T0OpuOSofIlDBly2nWHrzG281K06p6YdNx3MO9S7CkHWTJBx2WQ3prTw3UWjNu9zh2X9/NB3U/oI5PHdORhBAiySSmoG6utX6IY/nHVaA0MDxZUwnhBub9fpF5v180HSPVWxZ8manbzvKmXxH8Gz9jOo57iL4Hi9uALQ46rnIU1RY38/BM1p1dR98qfXnd93XTcYQQIkklZvu81+OPLwJLtdbh0g9WpAXfdfUzHSHVC/pTr+lPXpde0wAkxMKyTnDvAnReB3mtPzVw/dn1TD80nVdKvUL/KtbfVCmEEP8tMQX1RqXUSSAa6K+UygvEJG8sIczL5u31918knpr0mv4LWsN6f7i0E1rNgeLWnxq4O2w3Y3eNpY5PHcbWHSsXTUKIVOlv38G01qOAuoCf1joeiAReTe5gQpi2MTSMjaFhpmOkStfuR9Nj/j6yS6/p/7TtEziyAhqPgcrWnxp4KvwUQ4OGUiJHCaY0nIKXp1ykCiFSp799F1NKeQGdgece31n4FZiZzLmEMG7RnksAtKxS0HCS1MXRazqY6Dgbq6TX9P/bHwC/TYbqXaDB26bTuOxG5A36B/Yns1dmpjeZTtb00gZRCJF6Jea20Awc66inP37c+fFzPZMrlBDuYH73WqYjpDp/9Jq+cEd6Tf+Hs1vhh6FQqgm8NMXyvaYj4iIYEDiAyPhIAloEUCBzAdORhBAiWSWmoK6pta7yp8fblFKhyRVICHeRMb21e/66G+k17cSNI7CiG+QrD20DwOLLIuLt8QwNGsr5++f5tsm3lMlVxnQkIYRIdonZBWRTSv27+75SqiRgS75IQriHtQevsvbgVdMxUg3pNf0X7l+GRW+AdzbouAIyWPuOvdaaj3Z9xJ7re/ig7gfUK1TPdCQhhEgRiblDPRzYrpQ6DyigGNA9WVMJ4QaWBV8B4PVqUvy5avHeS9Jr+r9FhcOi1pAQDT1+gWzWX6s/M3Qm68+tp1+VftJrWgiRpvxtQa21DlRK+QJlcBTUJ7XWscmeTAjDFvWsbTpCqrD52A3GrDtKwzLSa/rf4qNhaXu4d9HRazpfOdOJXLbu7Dqmh07n1VKv0q9KP9NxhBAiRSWmy4c30B+oD2jgN6XUTK219KIWqZr0RXbd/kvhDFx6kEqFskuv6T/YbbC6J1zZC23mpYpe07vCdvHRro+o61OXD+t9KBdNQog0JzFLPhYAj4Cpjx+3BxYC1m+SKsQTrAxxLPlo41fEcBJrOnsrgrcCQvDJ7s3cbjXJlF56TaM1/DwCTv4ALT6FCtZfFnEq/BTDgoZRMkdJR69pD2tvqhRCiKeRmHe4Mv/V5WO7dPkQacGq/Y4NiVJQ/3M3H8bQdW4w6Tw8WNCjNrmzZDAdyT3s/BL2zYF6g6BOX9NpXHY94vq/e01/2+RbsqTPYjqSEEIYkZiC+qBSqo7Weg+AUqo28HvyxhLCvOV96pqOYEkPY+LpOjeY+1FxLO9Tl6K5M5mO5B4OLYXAj6BSW2j6kek0LnsQ+4C+W/sSHR/N/H/Nl17TQog0LTELGmsDu5RSF5VSF4HdwPNKqSNKqcNPe2CllKdS6qBS6ofHj3MppbYopc48/pjzaX+2EMKM2AQbvReEcPZWBDM716BioeymI7mHs4GwwR9KPA+vfgse1l5LHpMQg3+gP1ceXeHrxl9TOmdp05GEEMKoxNyhbpFMxx4MnACyPX48CgjUWk9SSo16/HhkMh1biL+1NPgyAO1rFTWcxBrsds2wFaHsOR/Ol29WoYFvXtOR3EPYIVjRBfKWgzcXQbr0phO5JMGewPAdwwm9Hcrk5ydTs0BN05GEEMK4xLTNu5TUB1VKFQZeAsYDwx4//SrQ8PG/A4AgpKAWBv1wOAyQgjoxtNZ88uMJfjx8nVH/Kiu9u/9w7yIsbgMZc0HHlY4BLhamtWb83vEEXQni3Vrv0rx4c9ORhBDCLZjadv8VMAL481iw/Frr6wBa6+tKqXxGkgnx2OKedUxHsIzvfjvP3N8v0P3Z4vR5rqTpOO4h8q5jcIs9Hjr9CNl8TCdy2YzQGaw6vYpelXrRoVwH03GEEMJtKK11yh5QqZeBF7XW/ZVSDYF3tNYvK6Xua61z/Onr7mmt/2cdtVKqN9AbIH/+/DWWLVuWUtH/Q0REBFmyyI52YT1J/drdFZbA7MOx1CrgSd8qGfCQHsR42GKpEjqGLBEXCK0yjofZrT+4ZeejnSwPX07tzLXpmLujkV7Tct4VViWv3dShUaNG+7XWfn/1ORN3qJ8FXlFKvQh4A9mUUouAm0opn8d3p32AW3/1zVrr2cBsAD8/P92wYcMUiv2fgoKCMHVskTIW7r4IQOe6xU3GSHJJ+drdeeYO87YEU6dkLgJ61CJDOs8k+bmWZot3TEF8dAbaLqB6uZamE7ks8HIgK4NW0qBQA75u/LWxXtNy3hVWJa/d1C/Ft5prrd/VWhfWWhcH2gHbtNadgA1A18df1hVYn9LZhPizrSdusfXEX17XCeDw1fv0WRhCqbxZmN3FT4ppALsd1g+As1vg5S8hFRTTB24eYOSOkVTIXYHJz0+WwS1CCPEX3Gl02SRghVLqLeAyMolRGBbQo5bpCG7r7K1HdJ0bTM7M6QnoUYts3lJkoTVsHg2Hl0OTD6BGN9OJXHb23ln8t/njk9mHb5t8SyYv6SkuhBB/xWhBrbUOwtHNA631XaCJyTxCiL937X40nb8PxtPDg0Vv1SZ/Nm/TkdzDzimwZzrU6Q/1h/3917u5G5E36Lu1Lxk8MzCj6QxyestoACGEcMad7lAL4Vbm7rwAQI/6JQwncR93ImLpPGcvEbEJrOhTl+J5MpuO5B72z4fAcVD5TWg+Hiy+MfNB7AP6be1HRHwE81vMp3BWaYMohBBPYu1xXUIko13n7rDr3B3TMdzGo5h4us0LJuxBNPO61aScj7V7KieZ4xvgh6Hg2zxVTEGMio9iQOAALj28xFeNvqJsrrKmIwkhhNuTO9RCODGnq0yA+0NMvI2eASGcvP6I77r64Vc8l+lI7uHCDlj9FhSuCW0CwNPaa8njbfEMCxrGkTtH+OL5L6jjI73YhRAiMax9K0UIkezibXb8lxwg+GI4X7StQqMyMnMJcIwUX9oBcpWC9ssgvbU37NnsNt7d+S6/h/3Oh3U/pGmxpqYjCSGEZUhBLYQTs3ecY/aOc6ZjGGW3a0auOszWE7cY90oFXq1ayHQk93DnrGMKYsac0HkNZLL2Hfs/Ror/cvEXhtUYRivfVqYjCSGEpciSDyGcOHDpvukIRmmt+fjH46w5eI1hzUqnugE3T+3hdVj4uuPfnddCtoJm8ySBqQensvL0SnpU7EH3it1NxxFCCMuRgloIJ2Z2rmE6glFTt51l3u8X6fFsCQY2fsZ0HPcQFQ6LWkH0Pei2EfJY//+XgGMBfHfkO1r7tmZI9SGm4wghhCVJQS2E+B8Buy4yZctpWlcvzPsvlUNZvA1ckoh56FjmcfccdFwJBauZTuSytWfWMjlkMs2LNWdMnTHyexZCiKckBbUQTkwPOgtA/4bWvwv5T6w/dI0PNxyjabn8fNq6Eh4eUmQRFwVL28GNw/DmIij5vOlELgu8HMjY3WOp61OXiQ0m4ukho+OFEOJpSUEthBPHwx6ajpDitp+8xdsrQqldIhfTOlQjnafsWyYhDlZ0gUu7oPUcKPMv04lctvf6Xob/OpyKeSryVaOvSO+Z3nQkIYSwNCmohXBiWofqpiOkqH0Xw+m7aD9lfbIyp6sf3l5yxxJbAqzpCWe3QMtvoNIbphO57MjtIwzaNohi2Yoxvcl0MnlZu92fEEK4A7n9JIQg9Mp9us/bR6EcGZnfvRZZva09oCRJ2O2wcRAcXw8vTIAaXU0nctmp8FP02dqHnN45mdl0JtkzZDcdSQghUgUpqIVw4pvAM3wTeMZ0jGR3POwhXeYGkzOzF4t71SZPlgymI5mnNWwaBYcWQ8N3oe4A04lcdv7+eXpv6U2mdJn4/oXvyZ85v+lIQgiRasiSDyGcOH87wnSEZHf21iM6f7+XjF6eLOlZB5/sGU1Hcg/bPoHgWVDXH54faTqNyy4/vEzPzT1RKOY0n0OhLDKgRwghkpIU1EI48VU767dFe5KLdyLp8N1elFIs6VWbIrlkLS0AO7+E3yZD9a7Q/BOweCu56xHX6bm5J/H2eOa+MJfi2YubjiSEEKmOFNRCpEHX7kfTcc5e4m12lvWuS8m8WUxHcg97Z8HWsVDxDXj5S8sX07ejbtNzc08i4iKY88IcfHP6mo4khBCpkhTUQjgxZfMpAIY1L2M4SdK6F2Onw3d7eBgTz9JedShTIKvpSO5h3xz4eQSUfRlenwkW78scHhNOr829uB19m9nNZlM+d3nTkYQQItWSgloIJ8IexJiOkOTuRsTy+b4YHsR7sLBnbSoWki4PAOyfDz++DaX/BW/MA09rdzl5EPuAPlv6cDXiKjOazqBqvqqmIwkhRKomBbUQTkxuU8V0hCR1LzKOTt8Hcydas6BnTaoXzWk6kns4uAg2Dgbf5tA2ANJZe8hJZHwk/bf25+z9s0xtPJWaBWqajiSEEKmeFNRCpAH3IuPoOGcv525HMKh6BuqUzG06knsIXQbr/aFUY2i7ENJZu2VgVHwUAwIHcOzuMb5o+AX1C9U3HUkIIdIE6UMthBOfbjrJp5tOmo7hsvtRcXT6fi9nb0fwXRc/KuaR62gADq+Edf2gRANotwS8vE0ncklUfBT+2/w5eOsgE+pPoEnRJqYjCSFEmiEFtRBO3I+K435UnOkYLvmjmD5zK4LZnWvwfOm8piO5h6NrYG1vKFoP2i8HL2v33/6jmN5/cz/j64/nxZIvmo4khBBpityqEsKJia0qm47gkgdR8XT6fi+nb0Qwu0sNGpbJZzqSezi+Hlb3hCK1ocNySG/t/tvRCdEM3Dbw38X0yyVfNh1JCCHSHLlDLUQq9OdielZnKab/7eSPsKoHFKoBHVdCBmv3345OiMY/0J+QmyF88uwnUkwLIYQhUlAL4cT4H48z/sfjpmP8Yw+i4+k8dy+nbjxiZufqNCorxTQApzbBiq7gUwU6rYIM1u6//d/FdMtSLU1HEkKINEuWfAjhREy83XSEf+xBdDydv9/LiesPmdmpBo3L5jcdyT2c+AFWdoP8FaDTGvC2dv/t6IRoBgYOZN+NfYyvP16KaSGEMEwKaiGc+Pi1iqYj/CP3o+LoOjeYE9cfMqNjDZqUk2IagGPrYPVb4FMVOq2GjDlMJ3LJH8V08I1gKaaFEMJNSEEtRCpwNyKWTt8Hc+5WBNM71qBpeSmmATiyCtb0hsJ+0HEVeGczncglf2xAlGJaCCHcixTUQjjx0cZjAHzYsoLhJE9262EMHebs5Up4FN919ZPWeH8IXeboM120LnRYYfkNiJHxkfgHOlrjfVJf1kwLIYQ7kYJaCAu7dj+ajt/t4dajWOZ3r0XdUjIBEYADC2HDQMfQlvbLIH1m04lc8jDuIf239ufonaNMajBJ+kwLIYSbkYJaCCfc/c705btRtP9uDw+j41n4Vm1qFMtpOpJ7CJkLPwx1jBNvt8TyQ1vux9yn95benLl/hi+e/4ImxWQCohBCuBspqIWwoLO3Iug4Zw+xCXaW9KpDpcLW7lqRZPbMhE0jwfcFaLvA8uPE70TfodfmXlx5dIVvGn1Dg8INTEcSQgjxF6SgFsKJMeuOAu7X7ePkjYd0mrMXgGW961C2gLU32iUJrWHH57B9PJR9Gd6YC+kymE7lkpuRN+m5uSc3o27ybZNvqe1T23QkIYQQTkhBLYQT3l7uN/foyNUHdJ67F+90nizuVZtSea290S5JaA2b34fd06BKe3hlGnha+9R2LeIaPX/pyb3Ye8xqNotq+aqZjiSEEOIJrP2uI0QyGv1SedMR/sP+S/foNjeY7Jm8WNKzDkVzZzIdyTy7DX4YAgcWQK3e0OJT8HC/C6F/4tLDS/Tc3JOo+CjmNJ9DxTzu9RcSIYQQ/0sKaiEs4NfTt+m7cD/5s2VgSa86FMxh7Y12SSIhDtb2gWNr4Lnh0Gg0KGU6lUtO3D1B3619AZj7wlzK5CpjOJEQQojEkIJaCCfeXXMYgImtKhvNsf7QNd5eEYpv/qwE9KhJvqzW3miXJOKiYGVXOLMZmn0Mzw4ynchl+27sY9C2QWRNn5VZzWZRInsJ05GEEEIkkhTUQjiRI1N60xGYu/MC4344Tp2SuZjdxY9s3l6mI5kX8xCWtoNLu+Dlr8Cvu+lELtt+eTvv/PoOhbMWZlazWRTIXMB0JCGEEP+AFNRCODGyRVljx9Za8/kvp5gedI4WFQrwVbuqeHt5GsvjNiLvwOI34MYRaD0HKr1hOpHL1p1dx9hdYymfuzzTm0wnh3cO05GEEEL8Q1JQC+FmEmx2Rq89yvKQK7SvVZRPXquIp4e11wYnifALsKgVPAyDNxdDmRamE7ks4FgAk0MmU9enLl81+opMXrLRVAghrEgKaiGceGdlKACT21RJsWPGxNsYuPQgW47fZFDjZxjarDTK4hvtkkTYIcedaXsCdNkARa3dk1lrzdcHvub7o9/TvFhzJjaYSHpP80uMhBBCPB0pqIVwomD2lN389yA6nl4BIey7FM5Hr1Sga73iKXp8t3VuGyzvDBlzQqc1kLe06UQusdltfLznY1afWU2b0m0YXXs0nh6ynEcIIaxMCmohnBjWPOValt16GEOXucGcux3B1+2q8UqVgil2bLcWuhzW94e8ZaHjKsjmYzqRS6ITohm1YxTbrmyjd+Xe+Ff1l79ACCFEKiAFtRCGXbgTSZe5e7kbEcfcbjVp4JvXdCTztIZdU2HLGCjeANotBu/splO5JDwmnIHbBnLk9hFG1RpFx3IdTUcSQgiRRKSgFsKJIcsOAvBVu+Qb+3zg8j16BYSggaW96lCliHR4wG6HzaNhz3So8Dq8PgvSZTCdyiWXH16m39Z+3Iy6yZcNv6RJsSamIwkhhEhCUlAL4UTJvFmS9edvCA3jnZWh+GT3Zm63mpRK5uNZQkIsrO3rmH5Ypz80H2/5UeKHbx/GP9AfjWZO8zlUzVfVdCQhhBBJTApqIZwY1MQ3WX6u1pqp284yZctpahXPxczONciVWTo8EBUOyzvBpd8d0w/rDbT8KPHtl7czYscI8mTMw4ymMyievbjpSEIIIZKBFNRCpKDYBBujVh9h7cFrtKpWiImtK5EhnXR44M4ZWNzG0WO69fepYmDLspPLmBg8kfK5yjOtyTRyZ8xtOpIQQohkIgW1EE74LzkAwLQO1ZPk54VHxtFnYQj7Lt7j7Wal8W/8jHR4ADgfBCu6gGd66PYDFKllOpFLbHYbk0Mms+jEIhoWbsinz30qA1uEECKVk4JaCCfKF8yWZD/r7K0Ieszfx42HMUxtX42W0hbPIWQe/Pg25C0D7ZdBzmKmE7kkIi6CETtG8Nu13+hUrhNv+71NOg85zQohRGonZ3ohnOjf8Jkk+Tm/n71Dv0X7SZ/Og2W961C9aM4k+bmWZrfB5jGw51t4phm8MRe8k+4CxoRrEdfwD/TnwoMLjKkzhrZl2pqOJP6vvTsPr6q69z/+XpnnQEICIYEEMIZR5hkVBK0zbcGxzkNbf7dW22rrxaqtra1Wq7deqqJWwWuLOLYidS5REGQeZEgChCkQCJCQeTg5Z/3+OAECEkjIsM9JPq/nyZMz7L3P98Tl4ZNv1l5bRKSNKFCLtKI3lu/i1//cQO+ESP5280h6xOlP/yVHTdoAACAASURBVFSXwjt3QM5HMPrH3pU8Av37o2htwVruWXgPLo+LFy58gTFJY5wuSURE2pB//ysm0op+/H+rAHjhxuFN3tfjsTzxURazvszl3PQu/PUHw4gJC27pEv3P4d0w91oo2AyX/RlG3uF0Rc02f9t8HlnyCEmRScycPJNesb2cLklERNqYArVIA4alntlFVipqarn3jbV8smk/N4zpyW+uGEBQoH+vpdwi8lZ5w3RtFfzgLTjLvy9u4rEeZq6ZyUvfvMTIbiN5+vyn6RSmC/OIiHRECtQiDfjheX2avM+uQxX8+PVVbN5XwsOX9+fW8WlayQNg3TyY/1OI6go3z4fEvk5X1CylNaU8uPhBFu5eyLT0aTw4+kGCA/UXCBGRjkqBWqSFfL55Pz+btxaAV24eyaS+iQ5X5APctfDpQ97LiKedC1fNhsguTlfVLLmHc7ln4T3sLt3NA6Me4Pq+1+uXJhGRDk6BWqQBd8xZAcDLN4885XZuj+WZT3OYuXArA7rH8MINw3XyIUD5QXjrFtixCEbfBRf9Dvy8i/v5zs+ZsXgGYUFhvHSRd6qHiIiIArVIA8b1OX0n9VBZNfe8sZbFWw9y9YgUHp06kLBgXfmQ/HXwxg+grAC++wIMuc7piprF7XHz17V/5aVvXmJg/ECemfQM3SK7OV2WiIj4CAVqkQbcNuHUqzWs2VXEf/19NQfLa3hi2iCuGdmzjSrzcWvnwgc/g4g4uO0jSG6ZK006pbi6mAcWPcDiPYv53lnf48ExDxIaGOp0WSIi4kMUqEWayFrL61/v5NEPNtE1Jox37xrHwORYp8tynqsKPvoVrJoNqRPgqlchyr/nkecU5XDvwnvJL8/noTEPcdXZV2m+tIiIfIsCtUgDbn5lOQBzbht19LGKmloefG8D763Zw6SMBJ65ZgidIkKcKtF3FO2AN2+G/LUw/l644CG/v1jLgtwF/Hbpb4kMjuTV77zKkMQhTpckIiI+yr//xRNpRVP6Hd9dzT1Qxl2vryanoJSfX3g2P5l0FgEB6laS8zG8+0OwFq6dC30vdbqiZql2V/P48sd5O+dthiYO5anznyIxwr877SIi0rraPFAbY3oArwHdAA/worX2L8aYOGAekAbsAK621ha1dX0iR9w4Nu3o7Y827OP+t9YRFGiYc+sozjs7wbnCfIW7FjL/CIuegm6D4OrXIK6301U1y66SXfzii1+QVZjFrQNv5e6hdxMc4N8rk4iISOtzokNdC/zCWrvaGBMNrDLGfArcAnxurX3cGPMA8ADwKwfqEzmq1u3hyU+ymfVFLoNTYnnuhuEkdwp3uiznFefBO3fArqUw9Ea49EkI9u+fy8c7PuaRJY8QaAKZecFMzu9xvtMliYiIn2jzQG2tzQfy626XGmM2A8nAVGBi3WZzgEwUqMVBV72wlOz9JZRU1nLDmJ48dHl/QoO0JB5ZC+Cf/w88tfD9l+Gcq5yuqFlq3DU8tfIp5mbN5ZyEc3jqvKdIikpyuiwREfEjxlrr3IsbkwZ8CQwEdllrO9V7rsha2/kk+/wQ+CFA165dh7/xxhttU+wJysrKiIqKcuS1pfWt2l/Li+urqfXAbQNDGJ/cfv7sf6Zj13hc9Nn2Kil7FlAa1YdN/e+nMsK/g2d+TT5zDs5hj2sPk6IncWXnKwkyOrXEV+lzV/yVxm77MGnSpFXW2hEne86xQG2MiQK+AB6z1r5rjDncmEBd34gRI+zKlStbu9STyszMZOLEiY68trSe4koXv52/kXdX72FgcgzPXD2E9K7RTpfVos5o7B7cCm/fCvvWw5j/gimPQJD/rsVsreWtnLd4csWThAeF8/sJv+e8lPOcLktOQ5+74q80dtsHY0yDgdqRVowxJhh4B/i7tfbduof3G2OSrLX5xpgkoMCJ2qTjWrzlIPe/vY6C0mp+Ojmduy84i+DAAKfLcpa1sPo1+Oi/vQH6unmQcbHTVTVLUVURjyx5hIW7FzKu+zh+P/73JEToJFMRETlzTqzyYYC/AZuttU/Xe+p94Gbg8brv/2rr2qRjqqxx8/iHm5mzdCd9EiJ5965xDO7RiWtmLQVg3o/GOlyhQ8oOwPyfQva/odd53kuIxyY7XVWzLMtfxoxFMyisLuS+EfdxY/8bCTAd/JcmERFpNic61OOBG4FvjDFr6x6bgTdIv2mMuR3YBfj3mU7iF9bsKuIXb64j92A5t43vxS8vziAs2Hvi4fThKQ5X56DsD+H9u6GqBL7zBxh9FwT4b/Csdlczc81M5mycQ2pMKjMnz6RffD+nyxIRkXbCiVU+FgMNXQ1jclvWIh1XTa2HZz/fwnOZW0mKDecfd45mXJ8ux21z1YgeDlXnoOoy+HgGrJ4DXQfBTe9D1/5OV9UsGw9t5MFFD7KteBvTz57O/SPuJyI4wumyRESkHdHp7NLhZO0r4efz1rEpv4Srhqfw0BX9iQn79ioeLrcHoOPMo969HN77ERRuh/H3wKQH/frEQ5fHxUvrX+LF9S8SHxbP81OeZ0LyBKfLEhGRdkiBWjoMt8fy8qJc/vxJDjHhQbx00wgu7N+1we1veHkZ0AHmUNdUwMLHYOlfITYFblkAaeOdrqpZthRt4cHFD7K5cDOX976cB0Y9QGxorNNliYhIO6VALR3CzkPl3PfWOlbsKOLiAd147HsDiY86dff12lEdYMrHjq/g/Z9AYS6MuB2m/AbCYpyu6oy5PW7mbJrDzDUziQ6J5pmJzzAldYrTZYmISDunQC3tmrWWfyzfxWMLNhMYYHjmmsF8d0gy3sVmTu17Q9vxSYnVZfD5b2H5i9ApFW6e713Jw4/lFOXwmyW/4ZuD3zC552QeGvMQ8eHxTpclIiIdgAK1tFv7iqv41Tvr+SLnABPO6sKfpp9D907hjd6/ssYNQHhI+7rceKeidfD83XB4N4z+MUx+GEIinS7rjFW7q5m1bhavbniV6JBoHj/3cS7tdWmjfmkSERFpCQrU0u54PJZ3Vufx+wWbqa518+jUAdwwOpWAgKYFrFteXQ60oznUlUXw6SMMWTcH4s+CWz+EVP9+byv2reDRpY+yo2QHV/a5kvtG3EfnsFNeYFVERKTFKVBLu7JqZxGPzt/Iurxihqd25snp59A7IeqMjnXDmNQWrs4h1sK6N+CTX0NlEbt6fJeeN70AwY3v1vuakpoSnl75NO9seYfkqGRmXTiLcd3HOV2WiIh0UArU0i7kF1fyxIdZ/HPtXrrGhPLMNYOZOji5yV3p+q4Y3L0FK3TIgWxY8AvYsQhSRsLl/yQ36xA9/TRMW2v5bNdn/GHZHyisKuSWAbdw1+C7tK60iIg4SoFa/FqVy81LX+byXOY23Nbyk0lncdfEPkSGNn9ol1S5AE66RrXPq6mARU/BV89650df/j8w7Gbv1Q6zMp2u7oxsL97OE8uf4Ku9X9Evrh8zJ89kQPwAp8sSERFRoBb/ZK3lww37eGzBZvYcruSSgd2YcWk/esS1XKfyzjkrAT+cQ53zCfz7Pji8EwZfBxf+DqISnK7qjFW4Kpi1fhavbXqNsMAwfjXyV1zb91qCAvTxJSIivkH/Ionf2bi3mEfnb2LZ9kL6dotm7p1jGNun5ZdHu3V8Wosfs1UV7fDOk948H7qc7fdL4Vlr+XjHxzy58kkKKgqY2mcq9w6/ly7hXU6/s4iISBtSoBa/caismqc+yWHeil3Ehgfz2PcGcu3IngQ2Y570qVw8MKlVjtviqkpg0Z/h6+cgIAgu+DWMuweCQpyu7IxtKdrCH5f/kRX7VtAvrh9/Pv/PDEkc4nRZIiIiJ6VALT6vptbDa0t38JfPt1BZ4+aWcb24Z3I6sRGtO7e5sLwGgLhIHw2mHjeseR3+8zsoP+Cd3jH5YYjx35MpS2tKeW7tc8zNmktUSBQPjXmIaenTCAxoX2uBi4hI+6JALT5tYXYBv/tgE7kHyjn/7AQeurw/ZyWe2TJ4TXXX66sAH51Dvf1L+GgG7P8GeoyG6+dB8nCnqzpjLreLednzmLV+FsXVxVx19lXcPfRuOoV1cro0ERGR01KgFp+07UAZv/9gEwuzD9C7SySv3jKSSX0T27SGO8/t3aav1yiHtsGnD0PWBxDbE6a/CgO+B356VUCP9fDxjo95dvWz5JXlMSZpDD8f/nP6xfdzujQREZFGU6AWn1Jc6eLZz7cwZ8kOwoMD+fVl/bhpbBohQQFtXsuU/l3b/DUbVH4IFj8Ny2ZBYAhc8BCM/S+/vjjL8vzlPL3qaTYe2khG5wxmTZnFuGRdnEVERPyPArX4BLfHMm/Fbp76JJuiihquHdmDX1yUQZeoUMdqKiitAiAxOsyxGqgqhqV/9X7VlMOQH8DkhyC6m3M1NdOWoi08s+oZFu1ZRLfIbjw24TEu63WZ5kmLiIjfUqAWx32de4jfzt/E5vwSRqXF8fAV/RmYHOt0Wdz9jzWAQ3Ooa8ph+Yuw+H+g6jD0nwoTZ0Bi37avpYXsK9/Hc2uf41/b/kVkUCQ/H/5zru93PaGBzv3SJCIi0hIUqMUR1lq+3HKQlxflsmjLQZI7hTPz+qFcNigJ4yPzge+a2KftX7S2GlbNhi+fgvICSL8IJj0I3f13ybiDlQd5ZcMrvJn9Jh7r4cZ+N3LnOXcSG+r8L00iIiItQYFa2lSVy837a/fy8uJccvaXkRAdyi8vzuC28b0IC/atP/lPzGjDkyDdtbDuH/DFn6B4N6ROgGv+D3qOabsaWlhhVSGzN8xmbtZcajw1XNH7Cu4achfJUclOlyYiItKiFKilTRwqq+bvy3bx2tIdHCyroW+3aP581WAuH5xEaJBvBekj9h6uBKB7p1Y88c9dC5v+CQv/AIXbvEvfXfks9J7ktyt3FFUV8dqm1/j75r9T7a7m0l6X8qNzfkRabJrTpYmIiLQKBWppVVsLyvjb4u28uzqP6loPEzMSuPPc3ozrE+8zUzsa8rN5a4FWmkPtqoS1/4Alz3ovGZ44AK6dCxmX+G2QLqgoYM7GObyV8xZVtVVcnHYxPx78Y3p38sHlB0VERFqQArW0OGstS3MP8fKi7fwnq4CQoACmDUvmtvG9SO8a7XR5jXb3Bektf9DKw7DiZVj2gvfqhsnD4aLfQ8ZlEND2SwO2hD1le3jlm1d4b+t7eKyHS3pdwh2D7qBPJwfmoIuIiDhAgVpaTE2thw/W7+XlRdvZlF9CfGQI905J54YxqY4uf3emJqR3abmDleyFr5+Dla9CTRmcNQXG3wtpE/y2I51TlMOcjXNYkLuAABPA1LOmctvA2+gR3cPp0kRERNqUArU02+GKGv6xfBdzluxgf0k16YlRPDFtEFOHJPvciYZNsetQBQA94yPO/CAHcmDJX2DdPLBuGDgNxt8D3Qa1UJVty1rL1/lfM2fjHL7a+xXhQeFc1/c6bh5wM90i/XdtbBERkeZQoJYztuNgOa9+tZ03V+ZR6XJzbnoXnph2DuefneDz86Mb4/631wFnOIc6byUsfgayFkBQKAy/Bcb9BDqntWiNbcXlcfHR9o+Ys3EO2UXZxIfF89OhP+XqjKu1/J2IiHR4CtTSJNZaVu4s4qUvc/l0836CAgxThyRz+4Re9EuKcbq8FvWzC89u2g611bDpfVj5N9i1FMI6wXn3wagfQVRC6xTZyg5VHuKdLe8wL3seBRUF9I7tzaPjHuXS3pfqgiwiIiJ1FKilUWrdHv69YR9/W5TLurxiOkUE818Tz+Kmsakkxjh4ae5WNKZ3fOM2LMz1XoxlzetQcQg694Lv/AGG3QSh/nMSZn0bDm5gbtZcPtz+IS6PizFJY3hk7CNMSJ5AgPHPkydFRERaiwK1nFJJlYt5y3cze8kO9hyupFeXSH733YFMH5ZCeIj/zo9ujG0HygDokxD17SfdtZDzEax8BbZ9DibQu+TdyNuh10S/XLGjxl3DJzs/Ye7muaw/uJ6IoAimpU/jur7Xaek7ERGRU1CglpPaXVjB7CU7mLdiN2XVtYzpHcdvrxzABX0TCQjw//nRjTHj3W+AE+ZQl+yF1a/BqjlQuheiu8PE//Z2o2O6O1Rp8+QW5/Juzru8v+19iqqLSItJ44FRDzC1z1SiQk7yy4SIiIgcR4FajrNmVxEvL9rOhxvyCTCGy89J4vYJvRmU0vFOPPvlxRneGx4PbM/0dqOz/u1draPPZLj0STj7Ygj0v/+Nqmqr+GTnJ7yT8w6rC1YTZIKY1HMS09OnM6b7GE3rEBERaQL/SwLS4qpcbhZmFfDy4u2s2llEdFgQd57Xm1vGpZEU24qX3fZxw+PdsG4u/OtV7zzpiHjvSh3Db4E4/5wCkVWYxZuFbzLjzRmUukrpGd2Tnw3/GVf2uZIu4S247raIiEgHokDdQe06VEFmTgGZ2QdYuu0QlS43PeLCeeSK/lw9ogeRoR10aFQe9i51t+FtsrduA+smI60HTJwB/a/0LoHnZ/aX72fB9gV8kPsBW4q2EEQQF/W6iOlnT2dE1xHtYolDERERJ3XQ1NTxVLncLN9eSGb2ATKzC8g9WA5AanwEV49IYVLfRM5NTyCwg8yPPk5NBeR8CBvehS2fgLsGOqfxcMQMiExg3m2Tna6wycpd5Xy28zPm585nef5yLJZzupzDjNEziNkTw2XnXeZ0iSIiIu2GAnU7truwgszsAhbW60KHBAUwtnc8N45NZWJGIr26RDpdpjNqq2Hbf+CbtyH7Q3CVQ3QSjLzTezXD5GHMyCt2usomcbldfJ3/NR/kfsDC3QuprK0kOSqZHw3+EZf1uoy02DQAMvdlOlqniIhIe6NA3Y4c14XOKSD3wPFd6IkZiYzpHd/ul7trkMcN27+EDe/A5vehqhjC4+Ccq2HQdOg5FgKO/WwG9+jkYLGN43K7WJq/lE92fMLC3QspqSkhJiSGy3tfzhV9rmBIwhBN6RAREWllCtR+7kgXOjP7AEvqdaHH9I7nhtGpTOrbgbvQ4F2hI28FbHgbNv4TygsgJBr6Xe7tRPeeCIHBJ911415vh3pAd99a4aTGXcOSvUv4dOenLNy1kFJXKVHBUUzqMYmL0i5iXPdxhASGOF2miIhIh6FA7Weqa491oRdmH+tC94xTF/oojxv2rIKsD7zzoot3Q1AYnP0dGDgd0i+E4NOvXvLo/E3ACetQO6TCVcHS/KV8tvMzMndnUuYqIzokmkk9J/GdtO8wJmmMQrSIiIhDFKj9wOm60BMzEujVJbJj/2m/eI/3ioVbP4fcTKg6DAFB3vWiL3jIexXDsJgmHfLhK/q3Tq2NtKdsD1/mfckXeV+wIn8FNZ4aYkJimJI6hYtSL2JM0hiCG+iui4iISNtRoPZB9bvQmdkFbKvXhb5qRAoTMxIY27tLx+5Cuyph5xLviYVbP4MDWd7Ho5O80zn6TPZO54iIO+OXaOupHm6Pm/UH1/PF7i/4Iu8Lth7eCkDP6J5c0/cazk85n2FdhxEcoBAtIiLiSxSofcTuwgoycw7wRXYBS7YdoqLG24Ue3SuOH6gLDdbCwRxvB3rrZ7DzK6itgsBQSB0HQ2/whujEftBCP6N1uw8DrXtyYmFVIcvyl/Fl3pcs3rOYw9WHCTJBDOs6jPtG3Mf5KecfXZ1DREREfJMCtUOqa92s2F7EwuyC47rQPeLCmT7c24Ue0zueiJAO/J+o8rB3+sa2z2Hrf6Akz/t4l7Nh+K1w1mRIHQ8hEa3y8n/492agZedQl7vKWbV/FV/nf82y/GXkFOUA0Cm0E+cmn8t5Pc5jXPdxxIQ0bXqKiIiIOKcDp7W2d6ou9PV1XejeHbkL7XHD3jXeLvS2z72rc1gPhMZC7/Pg/PuhzwXQqWeblPPo1IHNPkaNu4Z1B9axLH8Zy/KXseHgBmptLSEBIQxNHMpPh/6U0UmjGRA/gMCADjyFR0RExI8pULeiI13ozOwCMnMOsLWgDPB2oacNS2FS3w7ehfZ4oGg77FpadzLhQqgsAgwkD4Nz7/N2oZNHQGDb/4wyukU3eR+Xx0XWoSxW7F/BsvxlrN6/mip3FQEmgAHxA7hl4C2MThrNkIQhhAWFtULVIiIi0tY6aJJrPXlFFXUnEx5gybaD3i50YACje8dx3aieHbcLba13+bo9q71d6L2rYe86qK67GmFUN8i41NuB7j0JIuOdrRdYtbMQgOGpDZ/YWFhVyLqCdaw9sJa1BWvZeGgj1e5qAPrE9uH76d9ndNJoRnQboWkcIiIi7ZQCdTM11IVO6eztQk/MSGBsnw7YhS7dd0J4XgMVh7zPBQRDt4EwaBp0HwbJw1v0ZMKW8qePsoFjc6jdHjfbirextmAt6w6sY92Bdews2QlAUEAQ/eP6c3XG1QxJGMLQxKEkRCQ4VruIiIi0nQ6W8lpGXlEF/9nl4vU5K7/Vhb52ZA8mZiTSJ6EDdaHLD9UF5zXHAnRpvvc5E+gNyxmXQPeh3gDddQAEhTpbcyPMuDyNLUVbeX7t86w9sJb1B9ZT5vL+whQXFseQhCFMS5/G4ITB9I/vrykcIiIiHZQCdROVVrk4/8lM3B5LSueSjteFriqGvWuPD8+Hdx17Pj4dep13LDx3G9Rqq3C0FGste8v3klWYRXZh9tHve8v3AhBgAkjvlM6lvS5lSOIQhiQMISU6peP8wiQiIiKn1AESYMuKDgvm6asHU56XzXWXTWrfoaqmHPZ9c/zUjUNbjz3fKdU7XWPkHd7wnDS4yVcjbGsut4vc4lyyCrOOfmUXZlPqKgW84TktJo3BiYMZ7bqB1OhUrhkyjKiQKIcrFxEREV+lQH0Gpg5JJvPwlvYVpmurYf+GuvC81hueD2R5l60DiO7uXXlj8LXe8Nx9aLOuQtgWSmpKyC7MPtZ1Lspm6+Gt1HpqAQgPCie9czqX9LqEjLgM+sb1Jb1zOuFB4QBcM2sp2cDtoxSmRUREpGEK1B2FuxZK93qnZxzeXfd9Fxze6f1esgfqgiYRXbzhud8VdVM3hkJ0N2frb4DL42JP6R52luxkR8kOdpbsPHq7oKLg6HZdwruQEZfB+O7j6RvXl4y4DHpG9zzl2s9PTh/cFm9BRERE/JwCdXvhrvWG4iNBufiE0Fy8B6z7+H2ik7wXSekxyjt9I+kcb/c5NsWnVtzwWA8FFQXHheUjt/NK83DXe1+dQzuTGpPK2KSxpMWm0S+uHxlxGXQJ79Lk1+0Z79tzv0VERMQ3KFD7C7fr+MB89Gv3sQ7zcYHZ1AvMY2BQT+/tI1+xKT630kZxdfHRsLyj+Fho3lW6i8rayqPbhQWGkRqTSkbnDC5KvYi02DRSY1JJi0kjNjS2xepZvOUgABPSmx7GRUREpONQoPYVtTVQknfCdIx6X6V7j81nBsBATDJ06gGpY48Py516QkwKBIU49nZO5PK4OFhxkP0V+71f5fspqChgf8V+8svz2Vmyk8PVh49uH2gCSYlOITUmlVFJo0iL8Ybm1JhUEiMSCTABrV7z//5nC6BALSIiIqemQN1WaquhOO8k0zHqvkr2AvbY9iagLjD3hF7nQmyPEwJzss8E5gpXBfsrjgXkgooC9pXvO+7+ocpD2PrvDwgNDKVrRFe6RXbjwtQLj3aZU2NSSY5OJjgg2KF35PXMNUMcfX0RERHxDwrULeVoYN757ekYh3fVXejkxMCcUheYz68XlnscC8yBzgZKay2Hqw8fDcYndpYLKgrYX77/6JJz9cWExNA1siuJEYn0i+tHYkQiXSO897tGdqVrRFdiQmJ8eqWU7p3CnS5BRERE/IACdVO5XbDmdXrlLoaD/3es23zkyoBHmECITfae7Ndn0gnzl3tATPc2D8zWWiprKymtKaWkpuS478XVxd7OcsW+o0G5oKKAGk/NcccIMAF0CetCYkQiqTGpjOw28mhQ7hbZjcSIRBIjEo8uPefPMrO9q4RMzEh0uBIRERHxZQrUTWUC4cNf0sPj9naTY3tAn8nfnsMcnQSBLf/jdXlc3iBc7Q3CR0LxiQH5yPcTHzuyBvPJhASEHO0gD0oYRLeIbkfvH+kwdwnvQlBAxxg2z2duAxSoRURE5NQ6RjJqSQEBcO83fLlyIxMnTW7y7h7rodxV/q2geyQgNxSEj9yuv9rFyQQHBBMTEkN0SDQxITHEhMaQEpVCTKj3sSOPH32+3u3Y0FifnoLR1v73+qFOlyAiIiJ+QIH6DFSFd6LYXUbu4dyTd4erSyh1lR69XX+bMlcZnuNW6ziewRAVEnVc2E2LSTsuCEeHRBMTGnPcNkeeDw0MVShuIYnRYU6XICIiIn5AgbqJ3B43I/8+0ntnz8m3CQ8KJzo4+mhXOCEigT6d+jTcHa7XPY4KjmqTJeHk9D7btB+AKf27OlyJiIiI+DKfC9TGmIuBvwCBwMvW2scdLuk4gQGB3DfiPnbn7mb4wOHHBeQjoTgk0DeWs5PmeWlRLqBALSIiIqfmU4HaGBMI/BW4EMgDVhhj3rfWbnK2suPdPOBmMg9kMrHXRKdLkVb0/A3DnS5BRERE/ICvzS0YBWy11uZaa2uAN4CpDtckHVRcZAhxkfprg4iIiJyaT3WogWRgd737ecDo+hsYY34I/BCga9euZGZmtllx9ZWVlTn22tI2Vu7zLjE4opuv/W/SPBq74q80dsVfaey2f76WFE62PMVx16u21r4IvAgwYsQIO3HixDYo69syMzNx6rWlbTw/aykA91071uFKWpbGrvgrjV3xVxq77Z+vBeo8oEe9+ynAXodqkQ7upZtHOF2CiIiI+AFfm0O9Akg3xvQyxoQA1wLvO1yTdFAxYcHEhLXt5eFFRETE//hUh9paW2uM+QnwMd5l816x1m50uCzpoOav8/5x5IrB3R2uRERERHyZTwVqAGvtv4F/O12HyOtf7wQUqEVEROTUfC5Qi/iK2beOcroEERER8QMK1CINCA8JdLoEERER8QO+dlKiiM94b00e763Jc7oMERER8XHqzQW9ggAABj9JREFUUIs04I3l3msMfW9oisOViIiIiC9ToBZpwOt3jD79RiIiItLhKVCLNCA4UDOiRERE5PSUGEQa8NbK3by1crfTZYiIiIiPU6AWacDbq/J4e5VOShQREZFTM9Zap2s4Y8aYA8DOeg/FAsWN2LUx251umy7AwUa8lr9o7M/On167JY57psdo6n4au83j1PjV2NXYbS6N3ZY7hsZu22pvuaExx0y11iac9Blrbbv5Al5sqe1Otw2w0un368TPzp9euyWOe6bHaOp+GrvO/7f2pdfV2G3weY1dH39djd0Gn9fY9fHXbu4x29uUj/ktuF1jj9VeOPl+W+u1W+K4Z3qMpu6nsds8Tr1njV2N3ebS2G25Y2jstq32lhuadUy/nvLhJGPMSmvtCKfrEGkqjV3xVxq74q80dtu/9tahbksvOl2AyBnS2BV/pbEr/kpjt51Th1pEREREpBnUoRYRERERaQYFahERERGRZlCgFhERERFpBgVqEREREZFmUKBuBcaY7xpjXjLG/MsYc5HT9Yg0ljGmtzHmb8aYt52uReRUjDGRxpg5dZ+1P3C6HpGm0Gdt+6NAfQJjzCvGmAJjzIYTHr/YGJNtjNlqjHngVMew1v7TWnsncAtwTSuWK3JUC43dXGvt7a1bqcjJNXEMfx94u+6z9so2L1bkBE0Zv/qsbX8UqL9tNnBx/QeMMYHAX4FLgP7AdcaY/saYQcaYD074Sqy366/r9hNpC7NpubEr4oTZNHIMAynA7rrN3G1Yo0hDZtP48SvtTJDTBfgaa+2Xxpi0Ex4eBWy11uYCGGPeAKZaa/8IXH7iMYwxBngc+NBau7p1KxbxaomxK+KkpoxhIA9vqF6LmkPiA5o4fje1bXXS2vQh1DjJHOuEgPeDPPkU298NTAGmG2N+3JqFiZxGk8auMSbeGPMCMNQY89+tXZxIIzQ0ht8FphljngfmO1GYSCOcdPzqs7b9UYe6ccxJHmvwEpPW2meBZ1uvHJFGa+rYPQTol0DxJScdw9bacuDWti5GpIkaGr/6rG1n1KFunDygR737KcBeh2oRaQqNXfF3GsPizzR+OwgF6sZZAaQbY3oZY0KAa4H3Ha5JpDE0dsXfaQyLP9P47SAUqE9gjJkLLAUyjDF5xpjbrbW1wE+Aj4HNwJvW2o1O1ilyIo1d8Xcaw+LPNH47NmNtg9MpRURERETkNNShFhERERFpBgVqEREREZFmUKAWEREREWkGBWoRERERkWZQoBYRERERaQYFahERERGRZlCgFhFpB4wxjxpjpjTzGPcaY26qu/2EMWa9Mea1es/faIy5p979QcaY2c15TRGR9kDrUIuI+DljTKC11t2cfYwxQcBqYBgQCXxgrT3XGPN34HFgK/ABcLG11lVvv8+A26y1u1rgrYiI+CV1qEVEfJQxJs0Yk2WMmVPXLX7bGBNR99wOY8zDxpjFwFXGmNnGmOl1z002xqwxxnxjjHnFGBN6sn1OeLkLgNV1V3bzACHGGAOEAy7gfuDZ+mG6zny8l1MWEemwFKhFRHxbBvCitfYcoAT4f/Weq7LWTrDWvnHkAWNMGDAbuMZaOwgIAu461T51xgOrAKy1pcA7wBpgO1AMjLTW/usk9a0Ezm3G+xMR8XsK1CIivm23tfarutuvAxPqPTfvJNtnANuttTl19+cA551mH4Ak4MCRO9baP1lrh1hrfwH8DnjYGHOHMeZNY8yv6+1XAHRv/NsREWl/FKhFRHzbiSe61L9ffpLtzWmOd7J9ACqBsG8dzJihdTdzgJustVcDA40x6XWPh9XtKyLSYSlQi4j4tp7GmLF1t68DFp9m+ywgzRhzVt39G4EvGvE6m4GzTvL474CHgWAgsO4xDxBRd/tsYEMjji8i0m4pUIuI+LbNwM3GmPVAHPD8qTa21lYBtwJvGWO+wRt+X2jE63zI8VNDMMZ8F1hhrd1rrT0MLK07prXWrqvbbBKwoClvSESkvdGyeSIiPsoYk4Z3+bqBbfR67wG/tNZuaeT2oXi73xPqVgcREemQ1KEWEZEjHsB7cmJj9QQeUJgWkY5OHWoRERERkWZQh1pEREREpBkUqEVEREREmkGBWkRERESkGRSoRURERESaQYFaRERERKQZ/j8Qwlw5XdF1YwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = frm.plot(x='prior (%)', grid=True, figsize=(12, 7), \n",
    "              title='Pr(real Positive | test=Positive) @ Sensitivity = 99%',\n",
    "             logx=True)\n",
    "ax.set_ylabel('posterior (%)')\n",
    "ax.axvline(x=0.05, linestyle=':')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Q2. Given a test result is negative, what's the probability that the patient is really negative?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$ Pr(N|test=N) = \\frac{Pr(test=N|N) * Pr(N)}{Pr(test=N)} = \\frac{Pr(test=N|N) * Pr(N)}{Pr(test=N|N) * Pr(N) + Pr(test=N|P) * Pr(P)}$$\n",
    "\n",
    "Therefore,\n",
    "\n",
    "$$ Pr(N|test=N) = \\frac{Specificity * Pr(N)}{Specificity * Pr(N) + (1-Sensitivity)*Pr(P)} $$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "def posterior_neg(sen, spe, prior):\n",
    "    return spe * (1-prior) / (spe * (1-prior) + (1-sen) * prior)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.999900965311718"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "posterior_neg(.803, .995, 0.0005)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Specificity Study"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Following is quoted from the [preprint](https://www.medrxiv.org/content/10.1101/2020.04.14.20062463v1.full.pdf):\n",
    "\n",
    "\"The test kit used in this study (Premier Biotech, Minneapolis, MN) was tested in a Stanford laboratory\n",
    "prior to field deployment. Among 37 samples of known PCR-positive COVID-19 patients with positive\n",
    "IgG or IgM detected on a locally-developed ELISA test, 25 were kit-positive. A sample of 30 pre-COVID\n",
    "samples from hip surgery patients were also tested, and all 30 were negative. The manufacturer’s test\n",
    "characteristics relied on samples from clinically confirmed COVID-19 patients as positive gold standard\n",
    "and pre-COVID sera for negative gold standard. Among 75 samples of clinically confirmed COVID-19\n",
    "patients with positive IgG, 75 were kit-positive, and among 85 samples with positive IgM, 78 were kitpositive. Among 371 pre-COVID samples, 369 were negative. Our estimates of sensitivity based on the\n",
    "manufacturer’s and locally tested data were 91.8% (using the lower estimate based on IgM, 95 CI 83.8-\n",
    "96.6%) and 67.6% (95 CI 50.2-82.0%), respectively. Similarly, our estimates of specificity are 99.5% (95\n",
    "CI 98.1-99.9%) and 100% (95 CI 90.5-100%). A combination of both data sources provides us with a\n",
    "combined sensitivity of 80.3% (95 CI 72.1-87.0%) and a **specificity of 99.5% (95 CI 98.3-99.9%)**\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Now let's consider a similar toy problem:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let the random variable $X$ be the test result of a negative sample. $X=\\{0, 1\\}$, 1 means positive, 0 means negative, and $Pr(X=1)=p$.  So $X$ follows a Bernouli distribution $X \\sim B(p)$  \n",
    "\n",
    "In order to estimate $p$, we drew N=400 samples (i.i.d.), and found 2 of them are positive, ie, $X_i=1$ in two cases"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Q1: Let the sample mean be $p$'s estimator, ie,  $\\hat{p} = m = \\frac{1}{N}\\Sigma_i X_i$.  What is the expectation and variance of $\\hat{p}$ ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### A1:\n",
    "\n",
    "From the example, we can estimate $ \\hat{p} = m = 2/400 = 0.5\\% $.  Now let's check the expectation this estimate.\n",
    "\n",
    "$$E[\\hat{p}] = E\\big[\\frac{1}{N} \\sum_i X_i\\big] = \\frac{1}{N} \\sum_i E[X_i] = p$$\n",
    "\n",
    "So the estimator $\\hat{p}$ is unbiased. \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next we check the variance of the estimator:\n",
    "\n",
    "$$Var[\\hat{p}] = Var\\big[\\frac{1}{N} \\sum_i X_i\\big] $$\n",
    "\n",
    "Due to i.i.d, all the covariance terms are zero, so\n",
    "\n",
    "$$ Var[\\hat{p}] = \\frac{1}{N^2} \\sum_i Var[X_i] = \\frac{Var[X_i]}{N} = \\frac{p(1-p)}{N} $$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the sample variance is a consistent estimator of $Var[X_i]$, so we can calculate the sample variance of X\n",
    "$$s_X^2 = \\frac{\\sum_i (X_i-m)^2}{N-1} =\\frac{0.005^2 * 398 + 0.995^2*2}{399} = .004987$$\n",
    "\n",
    "Then the estimator s.d. of $\\hat{p}$ is\n",
    "$$s_p = \\sqrt{s_X^2/N} = 0.00353 $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the estimate of p is $\\hat{p} = 0.5\\%$ and $s_p = 0.353\\%$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "s.d. of p_hat 0.0035311006328336203\n"
     ]
    }
   ],
   "source": [
    "s2 = (.005**2 * 398 + .995**2*2) / 399\n",
    "print('s.d. of p_hat', np.sqrt(s2/400))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Q2. What's the 95% CI of $\\hat{p}$ ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Asuming $\\hat{p}$ follows a normal distribution $\\hat{p} \\sim N(m, s_p)$, then the 95% CI is $0.5\\% \\pm 1.96*0.35\\%$ \n",
    "\n",
    "Of course, the above goes out of bound, and we know $\\hat{p}$ doesn't follows a normal distribution.  \n",
    "\n",
    "A proper way is to use a bootstrap (resampling with replacement)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 127,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.011859999999999999 -0.0018599999999999997\n"
     ]
    }
   ],
   "source": [
    "print(.005 + 1.96 * .0035, .005 - 1.96 * .0035)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 195,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "index shape:  (1000, 5000) , index min: 0 , index max: 399\n"
     ]
    }
   ],
   "source": [
    "N = 400\n",
    "x = np.zeros(N)\n",
    "x[0] = 1   # should not matter as far as two of them are 1\n",
    "x[1] = 1\n",
    "\n",
    "batch_size = 5000\n",
    "resample = 1000\n",
    "idx = np.random.uniform(0, 400, size=(resample, batch_size)).astype(int)\n",
    "print('index shape: ', idx.shape, ', index min:', np.min(idx), ', index max:', np.max(idx))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 196,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0050146 0.0018 0.0088\n"
     ]
    }
   ],
   "source": [
    "p_samples = np.sum(x[idx], axis=1) / float(batch_size)\n",
    "\n",
    "print(np.mean(p_samples), np.min(p_samples), np.max(p_samples))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 197,
   "metadata": {},
   "outputs": [],
   "source": [
    "cdf = np.zeros(batch_size)\n",
    "for i in range(batch_size):\n",
    "    if i == 0 or cdf[i-1] != 1.0:\n",
    "        cdf[i] = np.sum(p_samples <= float(i) / batch_size) / float(len(p_samples))\n",
    "    else:\n",
    "        cdf[i] = 1.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 198,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x2eeac4f7128>"
      ]
     },
     "execution_count": 198,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAZrklEQVR4nO3de3CV933n8ffXQugC4iosA4IIKL5gx4kTxaRpM1Zp6mJvWietm9rtbh0nDmUndMa7f8TuzM5mdvqP2+7upB67ocTx0M5kwjJ2aEmHLltfTh2PYwJxfQETjAy2dQwCXUDSkTgS0vnuH5IcIXQ5l+c5z3OOPq8ZjTnP83vO+fLj6DM//87vdx5zd0REpPRdE3UBIiISDAW6iEiZUKCLiJQJBbqISJlQoIuIlAkFuohImZgX1QvX19d7U1NTXtf29/ezYMGCYAuao9SXwVA/BkP9OLuf//znne6+YqpzkQV6U1MTR44cyevaRCJBS0tLsAXNUerLYKgfg6F+nJ2ZvT/dOU25iIiUCQW6iEiZUKCLiJSJyObQp3L58mWSySTpdHrGdosXL+b48eNFqip/1dXVNDY2UllZGXUpIjIHzBroZvY08EXgvLvfMsV5A/4GuBsYAL7q7q/lU0wymaSuro6mpiZGn3ZqfX191NXV5fMSRePudHV1kUwmWbduXdTliMgckM2Uy25g6wzn7wI2jv1sA76bbzHpdJrly5fPGOalwsxYvnz5rP+3ISISlFlH6O7+kpk1zdDkHuAffPR7eF81syVmttLdz+ZTUDmE+bhy+rvE3aWhEbr6Bznd2U/6cqbor3/0/DCX3z5X9NctN3OlH5uW17KxIfhZhiDm0FcDbRMeJ8eOXRXoZraN0VE8DQ0NJBKJK84vXryYvr6+WV9wZGQkq3b5OnfuHI888givvfYaVVVVrF27lscee4yNGzd+1Gb79u1s3bqVL33pS7zyyis8/PDDVFZW8txzz1FTU/NRu3Q6fdXfM05SqVSs65vJ4LDzzMkhftY+QmrIGYn6q/1fy29fhUwyB/rx7nWVfOWG+YE/bxCBPtUwdMpfLXffBewCaG5u9skbCI4fP57V3HiYc+juzp133skDDzzAs88+C8Drr79+1WtWVlZSU1NDXV0d+/bt41vf+hYPPvjgVc9XXV3NbbfdFkqtQSjFjRx96ct8/+XT/L9j5zjePsDdt6ykqb6WNUtrWbOslsU1xf8Q+siRIzQ3Nxf9dcvNXOnH+oVVXLe4OvDnDSLQk8CaCY8bgTMBPG8kXnzxRSorK9m+fftHxz75yU/i7uzYsYMXXniBdevWMX6np6eeeoq9e/dy8OBBnnvuOX7wgx9EVfqc8cQLrfzdS6eoX1jF9x9oZsuNDVGXROfJCm5ZvTjqMkqe+rEwQQT6fmCHme0BNgM9+c6fT/Q/fnyMt8/0TnluZGSEioqKnJ9z06pFfPt3bp6xzdGjR/n0pz991fF9+/Zx4sQJ3nrrLc6dO8emTZv42te+xkMPPcTLL7/MF7/4Re69996ca5LctPek2f3Ke/zep1bzv7/yyajLEYmVbJYt/hBoAerNLAl8G6gEcPedwAFGlyy2Mrps8ep5hzLw0ksvcf/991NRUcGqVavYsmVL1CXNSX/z/Eky7vyXL1wfdSkisZPNKpf7ZznvwDcDq2jMTCPpMOfQb775Zp555pkpz2nVSrRaz6fYe6SN/7h5LWuW1UZdjkjsaOv/JFu2bGFwcJDvfe97Hx07fPgwS5cuZc+ePYyMjHD27FlefPHFCKucezIZ59Fn32TB/Ap2bNk4+wUic1Cstv7HgZmxb98+Hn74YR577DGqq6tpamriO9/5Dh0dHXz84x/n+uuv54477oi61DnlX4+f48j7F/ire29lRV1V1OWIxJICfQqrVq1i7969Vx1/4oknpmy/e/fukCuSxIkO6qrm8eXbVkddikhsacpFSsLLrR18dsNyKiv0lhWZjn47JPbe7+qnrfsSn99YH3UpIrGmQJfY++m7XQD82q8o0EVmErtAH9+BWQ7K6e8SpTeSF1lSW8n6et08WGQmsQr06upqurq6yiIIx78Pvbo6+O9rmGveaOvh46sXax+AyCxitcqlsbGRZDJJR0fHjO3S6XRJBOX4HYskf+nLI5w418f2G9dHXYpI7MUq0CsrK7O6u08ikYj1NxhKcN4+28tIxrm1cUnUpYjEXqymXEQme7PtIgCfUKCLzEqBLrF29ExvaN8dLVJuFOgSa0c/7OGW1YuiLkOkJCjQJbbSl0c4eT7FLat0wwORbCjQJbZOtPcxknGN0EWypECX2Do2dseqmzVCF8mKAl1i6+iZHhbXVNK4tCbqUkRKggJdYuvYhz3cvGqRdoiKZEmBLrF0eSTD8fY+3QFeJAcKdIml1vMphoYz3LxKH4iKZEuBLrF09MMeQB+IiuRCgS6xdOxML7XzK1inr8wVyZoCXWLp7TO93LRyERXX6ANRkWwp0CWW3uvqZ8MKjc5FcqFAl9i5NDTC+b5B1i6rjboUkZKiQJfYSV4YAGCNAl0kJwp0iZ0PukcDXSN0kdwo0CV2FOgi+VGgS+x80D3AgvkVLFswP+pSREqKAl1ip617gDXLavUdLiI5UqBL7HzQPaDpFpE8ZBXoZrbVzE6YWauZPTrF+cVm9mMze8PMjpnZg8GXKnOBuyvQRfI0a6CbWQXwJHAXsAm438w2TWr2TeBtd/8E0AL8LzPTBKjkrCM1SPpyhrXLFegiucpmhH470Orup9x9CNgD3DOpjQN1NjrpuRDoBoYDrVTmhLZurUEXyde8LNqsBtomPE4Cmye1eQLYD5wB6oA/dPfM5Ccys23ANoCGhgYSiUQeJUMqlcr7WrlS3PrylTOj44D2k0dJnC2dj3ji1o+lSv1YmGwCfaqlBj7p8W8DrwNbgA3Av5rZT9y994qL3HcBuwCam5u9paUl54IBEokE+V4rV4pbX77x3EnM3uHLv30H1ZUVUZeTtbj1Y6lSPxYmmyFQElgz4XEjoyPxiR4EfuSjWoHTwI3BlChzyQfdA1y3qLqkwlwkLrIJ9MPARjNbN/ZB532MTq9M9AHwmwBm1gDcAJwKslCZG8bXoItI7mYNdHcfBnYAB4HjwF53P2Zm281s+1izvwA+Z2ZvAc8Dj7h7Z1hFS/nSkkWR/GUzh467HwAOTDq2c8KfzwB3BluazDXpyyO096YV6CJ5Kp1lBFL2khcuAbBmWU3ElYiUJgW6xMZH34O+VCN0kXwo0CU22nvSAKxcohG6SD4U6BIbZ3vSmMG1dVVRlyJSkhToEhvtPWnqF1ZRWaG3pUg+9JsjsdHem2bl4uqoyxApWQp0iY32njTXLVKgi+RLgS6xcbbnkkboIgVQoEss9A8O05se5rrFWuEiki8FusRCe+/YkkWN0EXypkCXWBhfg96gOXSRvCnQJRbO9miELlIoBbrEQnvP6Pe4XKdAF8mbAl1i4WxPmqW1lbqxhUgBFOgSC+09aa1wESmQAl1i4WyPdomKFEqBLrFwrjet+XORAinQJXLpyyN09Q+xUksWRQqiQJfIne8dBKBBI3SRgijQJXJnx5Ysag5dpDAKdImctv2LBEOBLpEb3yWqZYsihVGgS+Tae9LUVc1jYdW8qEsRKWkKdInc2Z5LWrIoEgAFukRudJeoAl2kUAp0iVx7r249JxIEBbpEKpNxOlNDrKiriroUkZKnQJdIXRgYYiTjCnSRACjQJVKdqSEA6hcq0EUKpUCXSHWmRrf9K9BFCpdVoJvZVjM7YWatZvboNG1azOx1MztmZv8WbJlSrjr6RgN9Rd38iCsRKX2z7uQwswrgSeC3gCRw2Mz2u/vbE9osAf4W2OruH5jZtWEVLOVlfIS+YqFWuYgUKpsR+u1Aq7ufcvchYA9wz6Q2fwT8yN0/AHD388GWKeWqIzXI/IprWFSjXaIihcom0FcDbRMeJ8eOTXQ9sNTMEmb2czP7k6AKlPLW2TfE8oXzMbOoSxEpedkMi6b6TfMpnufTwG8CNcBPzexVd3/niicy2wZsA2hoaCCRSORcMEAqlcr7WrlS1H154oM01e4l/+8ZdT+WC/VjYbIJ9CSwZsLjRuDMFG063b0f6Dezl4BPAFcEurvvAnYBNDc3e0tLS15FJxIJ8r1WrhR1X/7VGz9hXX01LS2fiayGIETdj+VC/ViYbKZcDgMbzWydmc0H7gP2T2rzT8DnzWyemdUCm4HjwZYq5agzNUj9Qq1wEQnCrCN0dx82sx3AQaACeNrdj5nZ9rHzO939uJn9X+BNIAM85e5HwyxcSl8m43T1a9u/SFCyWlrg7geAA5OO7Zz0+K+Bvw6uNCl349v+talIJBjaKSqR0bZ/kWAp0CUy2vYvEiwFukTml9v+FegiQVCgS2R+ue1fgS4SBAW6REbb/kWCpUCXyHT0ja5B17Z/kWAo0CUynakh6jV/LhIYBbpEprNvUCtcRAKkQJfIaNu/SLAU6BIJbfsXCZ4CXSKhbf8iwVOgSyS07V8keAp0iYR2iYoET4EukdD3uIgET4EukdC2f5HgKdAlEtr2LxI8BbpEQtv+RYKnQJdIaNu/SPAU6BIJbfsXCZ4CXSLRkRrUB6IiAVOgS9FlMk53/xD1dfoeF5EgKdCl6LTtXyQcCnQpuvFt/9olKhIsBboU3fi2f43QRYKlQJei07Z/kXAo0KXotO1fJBwKdCm6jj5t+xcJgwJdiq4jpW3/ImFQoEvRadu/SDgU6FJ0nX3aJSoSBgW6FN3olIsCXSRoWQW6mW01sxNm1mpmj87Q7jNmNmJm9wZXopQTbfsXCc+sgW5mFcCTwF3AJuB+M9s0Tbu/BA4GXaSUj/Ft/5pyEQleNiP024FWdz/l7kPAHuCeKdr9GfAscD7A+qTMdIxvKtKHoiKBy2Yh8GqgbcLjJLB5YgMzWw18GdgCfGa6JzKzbcA2gIaGBhKJRI7ljkqlUnlfK1cqdl8e6xwBINl6nET3O0V73bDpPRkM9WNhsgn0qRYL+6TH3wEecfeRmdYWu/suYBdAc3Ozt7S0ZFnmlRKJBPleK1cqdl9e/PcP4cjrfOHXN/Mr1y4s2uuGTe/JYKgfC5NNoCeBNRMeNwJnJrVpBvaMhXk9cLeZDbv7PwZSpZSNj7b9a8pFJHDZBPphYKOZrQM+BO4D/mhiA3dfN/5nM9sN/LPCXKby0bb/am37FwnarL9V7j5sZjsYXb1SATzt7sfMbPvY+Z0h1yhlRNv+RcKT1TDJ3Q8AByYdmzLI3f2rhZcl5aozNaTpFpGQaKeoFFVHn3aJioRFgS5F1alt/yKhUaBL0Wjbv0i4FOhSNNr2LxIuBboUjbb9i4RLgS5F09k3BOjm0CJhUaBL0WiXqEi4FOhSNB19Y1MuGqGLhEKBLkXTmdK2f5EwKdClaDpSg6yoq9K2f5GQKNClaDpTQ9Qv1Bp0kbAo0KVotO1fJFwKdCma871prl2kQBcJiwJdiqJ/cJiu/iEal9ZGXYpI2VKgS1F8ePESAI1LayKuRKR8KdClKNq6BwBYs0wjdJGwKNClKJIXRkfoazTlIhIaBboURVv3ANWV12jZokiIFOhSFMkLl2hcWqtNRSIhUqBLUbRdGNAHoiIhU6BLUbR1D2j+XCRkCnQJXc+ly/Smh1mzTCN0kTAp0CV0yQujSxa1qUgkXAp0CZ2WLIoUhwJdQje+qUgfioqES4Euofuge4C6qnksqa2MuhSRsqZAl9Cd6uhn/bULtQZdJGQKdAndux0pNtQviLoMkbKnQJdQ9Q8Oc7YnzYZrF0ZdikjZU6BLqE539gOwXiN0kdBlFehmttXMTphZq5k9OsX5PzazN8d+XjGzTwRfqpSidztSABqhixTBrIFuZhXAk8BdwCbgfjPbNKnZaeAOd78V+AtgV9CFSml6t6Ofaww+tlxr0EXCls0I/Xag1d1PufsQsAe4Z2IDd3/F3S+MPXwVaAy2TClV73akWLOslqp5FVGXIlL2sgn01UDbhMfJsWPT+TrwL4UUJeXj3fMpNqzQdItIMczLos1Ui4d9yoZmv8FooP/6NOe3AdsAGhoaSCQS2VU5SSqVyvtauVKYfZlx593zA3ys+lLZ/3vpPRkM9WNhsgn0JLBmwuNG4MzkRmZ2K/AUcJe7d031RO6+i7H59ebmZm9pacm1XgASiQT5XitXCrMv27oHuHzwRVo+dRMtt68N5TXiQu/JYKgfC5PNlMthYKOZrTOz+cB9wP6JDcxsLfAj4D+5+zvBlymlqHVshYuWLIoUx6wjdHcfNrMdwEGgAnja3Y+Z2fax8zuB/w4sB/52bHv3sLs3h1e2lIKT5/oAuL6hLuJKROaGbKZccPcDwIFJx3ZO+PNDwEPBlial7kR7imvrqli6QDeGFikG7RSV0Lxzro8brtPoXKRYFOgSipGMc/J8n6ZbRIpIgS6haOseIH05ww0KdJGiUaBLKE6MfyCqKReRolGgSyjeaR8N9I36Ui6RolGgSyhOnOtjzbIaFlRltZBKRAKgQJdQvHOuT/PnIkWmQJfADQ1nONXRrxUuIkWmQJfA/aK9l+GMc9PKRVGXIjKnKNAlcIdOdQNw+7plEVciMrco0CVwh0530bS8loZF1VGXIjKnKNAlUJmM87PT3WxetzzqUkTmHAW6BOoX7X30pofZvF7TLSLFpkCXQB06PXpvk83rNUIXKTYFugTq0KluGpfWsHpJTdSliMw5CnQJjLvzs/e6tbpFJCIKdAlM6/kU3f1DfFYfiIpEQoEugXn19Oj6c30gKhINBboE5tCpLq5bVM3aZbVRlyIyJynQJRDuzqHTo/PnYzcKF5EiU6BLIN7rGqCjb1DTLSIRUqBLIA6dGlt/rg9ERSKjQJdAHDrdTf3C+WxYsSDqUkTmLAW6FMzdOXSqS/PnIhFToEvBXm+7yJmeNL+6oT7qUkTmNAW6FOzx50+ytLaSL9+2OupSROY0BboU5PW2i7x4ooOHPr+ehbohtEikFOhSkMefP8mS2koe+FxT1KWIzHkKdMmLu/P3r7zHC784zzc0OheJBf0WSs5605d59Nk3OfBWOy03rODBX2uKuiQRQYEuOXqj7SI7fvgaZy6m+fO7buQbn1/PNddoqaJIHGQ15WJmW83shJm1mtmjU5w3M3t87PybZvap4EuVqGQyTuv5Pp58sZV7d75CJgN7//RX+dM7NijMRWJk1hG6mVUATwK/BSSBw2a2393fntDsLmDj2M9m4Ltj/5UYyWScvsFh+geHSY399+SFEereH/3a28HhDL2XLtN7aZiO1CDvd/XzftcAx8/20pseBuALNzXwP//gVpbUzo/yryIiU8hmyuV2oNXdTwGY2R7gHmBioN8D/IO7O/CqmS0xs5Xufna6J227MMB//T+v51V0+7lB9p/L79q5JOPOxUuX6egbpKNvkK7+IUYyfnXDQz+d8vr6hVU0La/lP9y6ktvWLuVTa5eyYcUC7QYVialsAn010DbhcZKrR99TtVkNXBHoZrYN2AZQ3bCen/ziTK71ApDJZDh5Ib9r55oFlcbiKuPGxcbiFfNYON+oroCaeUbVPBhKp6mpqQagwozaytFrFlYa1fMMGAK6IdVN8u3Rf1i5WiqVIpFIRF1GyVM/FiabQJ9qODZ5mJdNG9x9F7ALoLm52Q9/++4sXv5qiUSClpaWvK6VK6kvg6F+DIb6sTDZfCiaBNZMeNwITB4eZ9NGRERClE2gHwY2mtk6M5sP3Afsn9RmP/AnY6tdPgv0zDR/LiIiwZt1ysXdh81sB3AQqACedvdjZrZ97PxO4ABwN9AKDAAPhleyiIhMJauNRe5+gNHQnnhs54Q/O/DNYEsTEZFc6LtcRETKhAJdRKRMKNBFRMqEAl1EpEzY6OeZEbywWQfw/tjDxUDPhNOTH08+Vg90hlrg1DUEfe1s7aY7n8vx2fo27L6Mcz9Ody6O/ThdXUFfF/Z7Ur/bM5/P5vjH3H3FlM/q7pH/ALtmejz5GHCk2DWFce1s7aY7n8vx2fo27L6Mcz9m22dx6MdC+jKX68J+T+p3O5h+nO4nLlMuP57l8XTHwlTI62V77Wztpjufy/Fs+jZMce7H6c7FsR8Lec1crgv7Panf7ZnP53r8CpFNuRTCzI64e3PUdZQD9WUw1I/BUD8WJi4j9FztirqAMqK+DIb6MRjqxwKU5AhdRESuVqojdBERmUSBLiJSJhToIiJloiwC3cwWmNnfm9n3zOyPo66nVJnZejP7vpk9E3Utpc7MvjT2fvwnM7sz6npKlZndZGY7zewZM/vPUdcTd7ENdDN72szOm9nRSce3mtkJM2s1s0fHDv8e8Iy7fwP43aIXG2O59KO7n3L3r0dTafzl2Jf/OPZ+/CrwhxGUG1s59uNxd98OfAXQcsZZxDbQgd3A1okHzKwCeBK4C9gE3G9mmxi95d34TapHilhjKdhN9v0oM9tN7n3538bOyy/tJod+NLPfBV4Gni9umaUntoHu7i8B3ZMO3w60jo0kh4A9wD2M3tO0caxNbP9OUcixH2UGufTl2O0Y/xL4F3d/rdi1xlmu70l33+/unwM0nTqLUgu/1fxyJA6jQb4a+BHw+2b2XaLZkl1qpuxHM1tuZjuB28zsz6MpreRM9578M+ALwL3jt2uUGU33nmwxs8fN7O+YdNc0uVpWt6CLEZvimLt7P7qPaS6m68cuQOGTm+n68nHg8WIXU8Km68cEkChuKaWr1EboSWDNhMeNwJmIaill6sfgqC+DoX4MQKkF+mFgo5mtM7P5wH3A/ohrKkXqx+CoL4OhfgxAbAPdzH4I/BS4wcySZvZ1dx8GdgAHgePAXnc/FmWdcad+DI76Mhjqx/Doy7lERMpEbEfoIiKSGwW6iEiZUKCLiJQJBbqISJlQoIuIlAkFuohImVCgi4iUCQW6iEiZUKCLiJSJ/w/n5bNfbMeU1wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "frm = pd.DataFrame({'Cdf': cdf})\n",
    "frm.plot(logx=True, grid=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 199,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0032"
      ]
     },
     "execution_count": 199,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.searchsorted(cdf, 0.025) / float(batch_size)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 200,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0072"
      ]
     },
     "execution_count": 200,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.searchsorted(cdf, 0.975) / float(batch_size)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
